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VOLUME  I 


TECHNICAL  DESCRIPTION 


INTRODUCTION 


The  Solid  Propulsion  Optimization  Code  (SPOC)  performs  detailed 
preliminary  designs  of  a  large  variety  of  solid  propellant  rocket  motors. 
Dimensions  of  the  propellant  grai:.,  nozzle,  and  pressure  vessel  are 
adjusted  by  the  code,  along  with  propellant  formulation  and  burn  rate,  to 
produce  a  motor  design  that  meets  performance  requirements  and  satisfies 
\  design  constraints  and  operating  limits,  and  that  has  been  optimized  with 
respect  to  a  performance  parameter  selected  by  the  user  from  a  menu. 

This  volume  of  the  User's  Manual  -  Volume  I  (Technical  Description)  - 
gives  the  basis  for  the  code  computations,  analytical  developments,  logic 
flow  charts  used  in  verification  checks,  and  error  messages.  Volume  II 
(User's  Guide)  contains  the  input  and  output  dictionaries  and  their  accom¬ 
panying  illustrations,  along  with  other  input  instructions  needed  to  execute 
the  code.  Volume  III  (Program  Description)  contains  the  subroutine  descriptions 
and  flow  charts,  cross -indices  of  common  statements,  subroutines  and  call 
statements. 

SPOC  was  prepared  for  use  by  a  motor  designer.  The  user/designer 
controls  the  direction  taken  by  the  search  through  the  inputs.  Information 
used  in  the  code  must  be  provided  by  the  designer,  but  no  more  is  required 
than  what  must  already  be  accumulated  in  order  to  prepare  a  detailed  pre¬ 
liminary  design - which  is  what  this  code  will  produce.  It  is  not  intended 

that  this  code  replace  final  detailed  stress,  thermal,  and  combustion  sta¬ 
bility  analyses;  it  will  monitor  certain  stress,  thermal  and  stability  param¬ 
eters  ;n  the  search  for  an  cptimized  design  so  that  the  final  arrangement  is 
more  likely  to  pass  detailed  analyses.  SPOC  will  do  no  more,  nor  will  it 
do  any  less,  than  a  good  designer  will  do;  but  the  code  will  do  it  much  faster, 
thus  enabling  the  designer  to  examine  more  approaches  and  more  combina¬ 
tions  than  previous';  possible. 

The  code  is  operational  on  an  IBM  4341  and  a  CDC  6600  computer,  and 
the  two  versions  are  almost  identical.  Double  precision  statements  required 
on  the  IBM  computer  have  been  deactivated  in  the  CDC  version,  but  they 
remain  in  the  code.  Arc  sin  and  Arc  cos  functions  are  ARSIN  and  ARCOS  in 
the  IBM  version,  but  are  ASIN  and  ACOS  in  the  CDC  version;  however,  the 
latter  can  be  kept  for  an  IBM  version  if  the  H-extended  compiler  is  employed. 

All  coding  is  in  FORTRAN  IV  language. 

There  are  115  subroutines.  No  external  devices  (other  than  a  printer) 
are  used.  There  is  no  overlaying  of  the  code  structure.  Core  storage 
requirements  aie  118K  decimal  words  on  the  Thiokol  IBM  4341  computer  and 
?.70K  octal  words  on  the  AFRPL  CDC  6600  computer. 
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Execution  time  depends  on  the  user  end  the  problem,  end  so  there  is 
no  single  representative  time.  About  8  seconds  ere  required  for  e  single 
pass  through  the  analysis  routine  for  a  minimum  option  short  burn  time 
motor  .  About  19  seconds  a?  e  needed  for  a  problem  that  includes  propellant 
formulation  adjustments,  cost  and  propellant  structural  analyses,  impulse 
efficiency  calculation  and  ballistic  simulation  at  two  temperatures.  Typically, 
300  to  500  calls  to  the  analysis  routine  are  needed  to  achieve  an  optimum 
solution. 

Five  AFR PL -supplied  sample  problems  have  been  successfully  solved 
with  SPOC.  These  problems  demonstrated  many  of  the  features  of  the  code. 

Another  aspect  of  the  demonstration  phase  of  this  project  was  a  briefing 
for  representives  of  propulsion  and  prime  contractors  held  on  16  June  1981. 
This  final  version  of  the  User's  Manual  (three  volumes)  supersedes  the  draft 
version  distributed  at  that  briefing. 

Revisions  made  to  the  manual  in  the  future  will  consist  of  added  or 
replacement  pages,  identified  with  revision  level  and  effective  date;  the 
original  page  numbering  system  will  be  preserved.  A  revision  cover  sheet 
will  serve  to  transmit  the  changed/added  pages  and  as  a  record  of  the  changes. 
Revisions  to  the  code  will  be  transmitted  in  the  same  manner  if  the  changes 
are  not  too  extensive;  a  listing  of  the  changed/added  coding  will  be  distributed. 
More  extensive  changes  will  require  electronic  transmittal  through  the 
AFR  PL  computer  center. 


MOTOR  AND  PROBLEM  DEFINITION 

SPOC  includes  models  fox  five  propellent  grain  configurations,  three 
forward  closure  and  two  aft  closure  arrangements,  and  six  noszle  configura¬ 
tions.  Any  combination  of  grai..,  closure,  and  nozzle  may  be  selected  except 
that  «  Type  4  grain  (conocyl)  may  be  used  only  with  the  Type  1  forward  closure 
(ellipsoidal). 

The  analyses  performed  hy  SPOC  are 

Thermochemistry  Pressure  veosel  structural 

Ballistic  Nozzle  thermal  and  structural 

Propellant  structural  Trajectory 

Weight  Combustion  stability 

Cost  Impulse  efficiency 

Flexibility  has  been  provided  for  the  user  /designer  so  that  the  code 
may  be  tailored  to  enable  varied  problems  to  be  solved.  These  choices  are 
described  in  the  following  paragraphs. 

Propellant  Grain:  choose  one  from  those  illustrated  in  Figure  1. 

Type  1:  Star 

Type  2:  Double  web  wagon  wheel 

Type  3:  Finocyl  (slots  in  forward  end) 

Type  4:  Conocyl 

Type  5:  Cylindrioally  perforated  (CP) 

Nozzle;  choose  one  from  those  illustrated  in  Figure 

Type  1;  Thin  shell,  composite  structure  as  the  insulating 
ablative  and  support  structure. 

Type  2;  Thin  shell  support  structure  with  insert  and  abla¬ 
tive  insulator. 

Type  3;  One-piece  ablative;  supersonic  blast  tube;  con¬ 
stant  diameter  support  structure. 

Type  4;  One-piece  ablative;  supersonic  blast  tube;  reduced 
diameter  aft  section. 

Type  5;  Subsonic  blast  tube;  without  expansion  cone. 

Type  6:  Subsonic  blast  tuba;  with  expansion  cone. 
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Forward  Closure;  choose  one  from  those  illustrated  in  Figure  3. 

Type  1:  Ellipsoidal 

Type  2*  Flat  plate  with  closure  secured  with  retaining  ring 

Type  3:  Flat  plate  with  closure  integral  with  case 

Aft  Closure:  choose  one  from  those  illustrated  in  Figure  4. 

Type  1:  Ellipsoidal 

Type  2;  None  (aft  closure  formed  by  nossle  entrance 

section) 

Other  choices  that  must  be  made  to  define  the  problem  are: 

1.  A  propellant  formulation  may  be  input  and  adjusted  as 
part  of  the  optimisation  (FORMAD=T),  in  which  case  the  thermochemistry 
routines  are  entered  every  time  the  design  is  evaluated  (except  for  some 
internal  by-passes  to  reduce  execution  time).  Another  option  is  to  input  a 
formulation  but  not  adjust  It  (FORMIN=T),  in  which  ease  the  thermochemis¬ 
try  routines  are  entered  only  for  the  first  evaluation  in  order  to  obtain 
basic  propellant  charac  teristics  for  the  ballistic  simulation.  The  third 
option  is  for  the  user  to  input  the  appropriate  ballistic  parameters  rather 
than  having  the  thermochemistry  routines  calculate  them  from  a  fcrmula- 
tion  (PROPIN=T).  The  proper  combination  of  these  three  inputs  is  shown 


below  (all  default  to  F). 

MODE 

FORMAD 

FORMIN 

PROPIN 

(1)  Formulation  input  and  adjusted 

during  optimization 

T 

F 

F 

(2)  Formulation  input,  but  not 

adjusted 

F 

T 

F 

(3)  User  supplies  required 

propellant  characteristics 

F 

F 

T 

2.  Impulse  efficiency  may  be  input  by  the  user,  calculated 
internally  with  the  AFRPL  SPP  "empirical  model",  or  calculated  with  a 
user-supplied  model  which  must  be  installed  in  subroutine  USEREF. 
EFMDL=T  is  the  flag  to  show  a  user  model  has  been  supplied.  SPPETA-T 
is  the  flag  to  specify  the  SPP  model. 

3.  Propellant  burn  rate  is  calculated  internally  with  the 
Vielle  model,  or  with  a  aser- supplied  model  which  he  must  install  in  sub¬ 
routine  USERRB.  RBMDL=T  is  the  flag  to  show  a  user  model  has  been 
supplied. 
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Flat  Plate  -  Type  2  (Retaining  Ring) 

-  Type  3  (Integral  with  Case) 


Figure  3.  Forward  Closure  Configurations  Available  in  SPOC 


£ 


May  be 
Inhibit  _ 


t\\\\ 


Ellipsoidal  -  Type  1 


Figure  4.  Aft  Closure  Configurations  Available  in  SPOC 
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4.  The  propellant  face  on  the  forward  and  of  a  grain  with  a 
Type  2  or  Type  3  forward  closure  and  on  the  aft  face  of  a  grain  with  either 
a  Type  1  or  Type  2  aft  closure  may  be  inhibited  through  uae  of  FWD1NH-T 
or  AFTINH=T,  respectively. 

5.  Belli stic  simulations  will  be  performed  at  both  the  low 
temperature  and  high  temperature  conditions  if  different  values  are  input 
for  THl  and  TLO.  Propellant  structural  analysis  is  performed  at  a  differ¬ 
ent  temperature  (TPROP)  than  is  the  low  temperature  ballistic  simulation. 
Pressure  vessel  structural  analysis  is  performed  at  the  high  temperature 
condition.  If  THI  is  input  equal  to  TLO,  only  one  ballistic  simulation  is 
performed;  propellant  structural  analysis  is  still  performed  at  TPROP,  and 
pressure  vessel  structural  analysis  is  performed  with  the  result*  of  the 
single-temperature  ballistic  simulation  (i.  e. ,  pressure  not  adjusted  to  some 
high  temperature  condition). 

6.  The  optimisation  routine  will  adjust  user- specified  param¬ 
eters  in  order  to  meet  all  performance  requirements  and  satisfy  all  design 
constraints.  In  addition,  the  user  may  specify  another  parameter  to  be 
optimized  by  setting  ICHOZE  to  one  of  the  following. 

0;  None  (default  value) 

1:  Minimise  cost 

2 :  Minimise  total  motor  weight 

3:  Maximise  total  impulse 

4:  Maximise  total  impulse -to -total  weight  ratio 

5:  Maximise  burnout  velocity 

7.  There  are  36  parameters  (not  all  on  one  problem)  whose 
values  can  be  adjusted  by  the  optimisation  routine  PATSH  to  achieve  an 
optimum  design  (Table  1).  Each  of  these  must  be  specified  by  the  user  us 
"T"  (maintain  at  input  value)  or  "F"  (do  not  maintain  at  input  value,  but 
adjust  during  pattern  search).  Default  value  is  T  (do  not  adjust). 

8.  A  trajectory  simulation  (point  mass,  flat  earth,  ballistic 
trajectory)  will  be  performed  if  specified  by  the  user  (FTRAJ«T).  If 
ballistic  simulations  are  performed  at  two  temperatures  (TLO  and  THI), 
then  trajectory  simulations  are  performed  with  esch  of  the  resultant  thrust¬ 
time  histories.  In  addition,  the  user  must  select  a  trajectory  termination 
option. 

9.  Motor  cost  will  be  calculated  with  the  Tri-Services  cost 
model  or  with  a  user-supplied  model.  FCOST=T  is  the  flag  to  specify  the 
Tri-Services  model;  CSTMDL  is  the  flag  to  show  a  user-supplied  model 
has  been  provided. 

10.  Either  a  contoured  or  conical  nozzle  expansion  section  may 
be  specified  (CONTUR=T  or  CONTUR-F,  respectively).  If  a  conical  exit 
section  is  selected,  the  initial  half -angle  of  the  expansion  section  (ALFA) 
must  be  input  equal  to  the  exit  half -angle  (ALFAEX). 
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ADJUSTABLE  VARIABLES  AVAILABLE  IN  SPOC  (contd.) 
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1  and  Z) 

R2A5  GRAIN5  FR2A5  GRAIN5  Propcilar.t  port  radius  at  Plane  5  (Grain  S 


11.  Several  analyses  ere  by-pee«ed  completely  unless  the 
user  specifies  otherwise. 


(*)  Propellent  structural  analysis  (PSTRUC*T) 

(b)  Combustion  stability  (FSTABsT) 

(c)  Trajectory  simulation  (FTRAJ=T) 

(d)  SPP  impulse  efficiency  (SPPETA=T) 

(e)  Thermochemistry  (FORMAD=T  or  FORMIN=T) 

(f)  Cost  (FCOST*T) 


,  ,  .  1Z-  T,he  us«r  m*y  provide  model,  for  certain  p»r.mel.r.  that 

londed!  f  **  "  I"*:  A  fl*«  “  •«  »•  •‘■ow  .  user  mod.l  h>.  been 
®  into  a  specified  subroutine  (T=  model  has  been  supplied). 


Load  in 


^lag 

Subroutine 

Parameter  to  Be  Supplied 

RBMDL 

USERRB 

Propellant  burn  rate,  RATE  (in/sec) 

SEMDL 

USERSE 

Propellant  nominal  strain  endurance.  SENOM 
(in/in) 

EOMMDL 

USERRH 

Propellant  rheological  property  to  be  defined 
by  user.  EOM  (units  by  user) 

CSTMDL 

USERCS 

Motor  cost,  COST  ($  or  $/unit) 

EFMDL 

USE REF 

Impulse  efficiency,  ETAISP  (%  x  0.  01) 

+ 

RSPNSE 

Combustion  response 

♦IRSPNS  -  £  in  namelist  STABIN 


13.  If  a  combustion  stability  analysis  (Reference  7)  is  desired, 
the  user  can  select  one  of  five  combustion  response  models  (one  is  user- 
supplied)  and  cm  specify  at  how  many  modes  stability  margin  is  to  be 
calculated. 


COMPUTER  CODE  ARRANGEMENT 


The  computer  code  has  an  overall  organisation  similar  to  that  shown 
in  Figure  5.  There  are  three  major  subprograms  (MAIN,  COMP,  and 
PATSH),  whose  functions  are  listed  on  Figure  S.  Flow  through  the  program 
will  then  be  that  of  Figure  6.  MAIN  first  reads  and  initialises  various  pa> 
rameters  and  prints  a  narrative  summary  of  the  problem  as  defined  by  the 
user  in  subroutine  CHECKIN.  A  call  is  then  made  to  COMP  for  the  first 
time  in  order  to  calculate  performance  of  the  motor  with  user-supplied 
initial  values.  Initial  values  of  penalties  are  also  calculated  (in  COMP)  and 
all  output  is  printed,  after  which  the  PRINT  flag  is  turned  off.  MAIN  then 
calls  PATSH  to  adjust  specified  parameters  in  order  to  minimise  the  payoff 
parameter  and  penalties  (described  in  detail  in  the  PATSH  module  description 
in  another  report  section).  Each  time  PATSH  adjust  one  or  more  of  the 
specified  parameters,  COMP  is  called  to  calculate  motor  performance,  payoff 
and  associated  penalties.  PATSH  builds  a  pattern  and  makes  adjustments  to 
minimise  the  OBJ  function.  When  there  is  no  further  decrease  in  the  payoff 
and  penalties,  the  PRINT  flag  is  turned  on,  COMP  calculates  the  performance 
with  the  last  set  of  adjusted  parameters,  and  results  are  printed. 

The  executive  subroutine  COMP  sets  up  the  user-  or  PATSH- supplied 
inputs  for  the  various  analyses  and  simulations  and  passes  the  results  of 
early  analyses  to  later  calculations  when  they  are  needed  (Figure  7).  For 
the  first  pass  through  COMP,  where  all  analysis  inputs  are  furnished  by  the 
user,  the  inputs  are  read  in  the  specific  subroutine  to  which  the  data  applies. 

On  all  subsequent  passes  through  COMP,  the  input  data  are  either  constant 
at  the  user-input  value  or  are  updated  according  to  the  PATSH  adjustments. 
Write  commands  are  also  given  within  the  individual  subroutines. 

The  first  call  by  COMP  is  to  one  of  the  grain  dimension  verification  and 
setup  subroutines  (SETUP1  for  Crain  Type  1,  SETUP2  for  Crain  Type  2/  etc.); 
these  subroutines  verify  the  geometric  validity  of  the  incoming  dimension  set 
and  calculates  other  dimensions  needed  by  the  ballistic  simulation  module. 
Subroutine  NOZINP  is  called  to  perform  the  same  function  for  the  nossle. 

If  the  problem  involves  a  propellant  formulation,  subroutine  TCHEM  is 
called  next  to  perform  thermochemical  analyses;  results  of  the  calculations 
are  used  in  IMPEFF  (impulse  efficiency),  SEC2SB  (ballistic  simulation). 

NOZL  (noazle  thermal  and  structural  analysis),  and  E488M2  (combustion 
stability).  Subroutine  IMPEFF  is  called  next'to  furnish  a  value  for  impulse 
efficiency,  if  specified  by  the  user. 
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MAIN 


MAIN: 

COMP: 

PATSH: 

CHEKEN: 


ontMiN 

COMP 

PATSH 

- - — - 

COMP 

COMP 


Reads  control  inputs;  initialises  some  parameters;  controls 
printout;  calls  search  routine 

Executive  subroutine  passes  information  between  subroutines; 
calculates  some  penalties  and  overall  objective  function  (OBJ); 
provides  printout 

Adjusts  specified  parameters;  evaluates  changes 
in  objective  function  (OBJ) 

Checks  compatibility  of  problem  definition.  Prints  narrative 
description  of  problem. 


Figure  5.  Overall  Code  Organisation 
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START 


READ  INPUTS 


RUN  INITIAL  ESTIMATED 
DESIGN 

PRINT  RESULTS 


SUBROUTINE 

COMP 


TURN  OFF  PRINT  FLAG 


PERFORM  MINIMIZATION 


SUBROUTINE 

PATSII 

TURN  ON  PRINT  FLAG 


RUN  FINAL  DESIGN 
PRINT  RESULTS 


STOP 


SUBROUTINE 

COMP 


SUBROUTINE 

COMP 

Figure  6.  Overall  Code  Flowchart 


Figure  7.  Organization  of  Executive  Subroutine  COMP. 


pjiww??? 


Subroutines  SEClSB  i.nd  SEC2SB  make  up  the  bailie  tic  simulation 
module.  The  first  time  they  are  called,  the  input  ballistic  parameters 
have  been  set  up  (in  COMP)  for  a  grain  conditioned  to  high  temperature 
conditions.  When  the  ballistic  simulation  is  completed,  subroutine 
HITEMP  uses  the  results  to  calculate  certain  performance  parameters  and 
operating  conditions  associated  with  high  temperature  motor  operation  (e.  g. . 
design  pressures,  minimum  burn  time,  etc).  The  predicted  values  are 
compared  with  user  input  limits  and  appropriate  penalties  are  calculated. 

Next  COMP  sets  up  ballistic  parameters  for  a  simulation  with  the  grain 
conditioned  to  low  temperature,  and  then  SEClSB  and  SEC2SB  are  called 
again.  Results  of  the  low  temperature  simulation  are  analysed  in  sub¬ 
routine  LOTEMP  for  performance  parameters  and  operating  conditions 
associated  with  low  temperature  motor  operation.  If  the  user  wants  Qnly 
to  study  a  problem  via  ballistic  simulation  at  a  single  temperature,  the 
second  simulation  is  skipped  and  results  of  the  first  are  analysed  (by 
making  THI  =  TLO)  in  subroutine  ONETMP  (that  combines  the  calculations 
of  HITEMP  and  LOTEMP). 

Once  the  results  of  the  ballistic  simulations(s)  are  available,  nozzle 
thermal  and  structural  analyses  are  performed  in  subroutine  NOZL. 
pressure  vessel  structural  analyses  are  performed  in  subroutine  CASEAN,  and 
(if  specified  by  the  user)  prope.lant  structural  analyses  are  performed  in 
subroutine  PROPST.  The  user  may  also  command  a  trajectory  simulation. 
Subroutine  TRAJIN  acts  as  a  mini- executive  subroutine  to  control  the  tra¬ 
jectory  simulations  for  a  one-or  two-temperature  problem.  Motor  cost  is 
calculated  in  subroutine  COST,  and  combustion  stability  characteristics  are 
determined  in  subroutine  E438M2,  if  specified  by  the  user. 
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OPTIMIZATION  PROCESS 


SPOC  combines  computer  models  for  solid  rocket  motor  per¬ 
formance  prediction  and  design  analyses  with  a  numerical  parameter 
optimization  technique.  As  stated  in  Reference  1,  this  combination 
requires  an  understanding  of  both  areas.  The  following  discussion  was 
taken  from  Reference  1  because  approaches  taken  in  the  TACMOP  and 
SPOC  codes  are  very  similar,  even  though  the  codes  have  different  end 
objectives. 


In  order  to  eliminate  misinterpretation,  several  terms  used 
throughout  the  remainder  of  the  discussion  are  defined  below; 

Performance  requirement  -  A  measure  of  acceptable  system 
operation  in  accomplishing  its  intended  purpose.  For  solid- propellant 
rocket  performance  requirements  typically  include  such  items 

as  range,  ilocity,  or  payload  delivered  to  a  specified  end  condition. 

In  SPOC,  performance  requirements  are  expressed  as  total  impulse, 
impulse -to- weight  ratio,  etc. ,  as  well  as  the  ultimate  end-item  require¬ 
ments  listed  above;  however,  the  trajectory  simulation  in  SPOC  is  not 
intended  for  complex  maneuvering  trajectories,  and  so  SPOC  should  be 
used  in  conjunction  with  more  sophisticated  trajectory  simulations. 

Design  parameter  -  A  length,  angle,  or  material  property  used  in 
describing  a  particular  design,  such  as  propellant  grain  length,  case 
diameter,  nozzle  half  angle,  or  propellant  burning  rate. 

Design  constraint  -  A  limit  imposed  directly  or  indirectly  on  the 
allowable  values  of  a  design  parameter,  such  as  maximum  length,  maxi¬ 
mum  nozzle  divergence  angle,  maximum  propellant  web  fraction,  or 
minimum  port- to- throat  area  ratio. 

Operating  limit  -  A  maximum  or  minimum  acceptable  level  for  a 
condition  produced  by  motor  operation,  such  as  maximum  acceleration, 
minimum  pressure,  or  maximum  velocity. 

Payoff  -  The  quantity  selected  ss  the  maximized  or  minimized 
variable  during  the  optimization  process,  such  as  maximum  range.  In 
SPOC,  corresponding  payoffs  are  total  impulse,  motor  weight,  cost,  etc. 
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Penalty  function  -  A  function  corresponding  to  a  particular  per¬ 
formance  requirement,  design  constraint,  or  operating  limit,  having 
zero  value  when  the  requirement,  constraint,  or  limit  is  satisfied  by 
the  design  being  evaluated,  and  having  a  non-zero  value  proportional 
to  the  amount  of  violation  of  the  particular  requirement  when  it  is  not 
satisfied. 

Objective  function  -  A  single-valued  function  fcr  a  particular 
design  representing  both  die  payoff  value  and  any  non- zero  penalty 
function  values  associated  with  that  design. 

The  design  problem  consists  of  finding  a  set  of  design  parameter 
values  that  produce  a  system  with  maximum  (or  minimum)  payoff,  subject 
to  meeting  all  performance  requirements,  design  constraints,  and  operat¬ 
ing  limits  (i.  e.  ,  all  penalties  non- zero). 

Parameter  Optimization  Scheme 

The  optimization  routine  used  in  SPOC  is  the  PATSH  (Pattern 
Search)  subroutine  developed  by  D.  E.  Whitney  at  the  Massachusetts 
Institute  of  Technology  (Reference  2).  This  subroutine  performs  an  uncon¬ 
strained  non-linear  optimization  with  the  direct  pattern  search  algorithm 
of  Hooke  and  Jeeves  (Reference  3)  in  a  highly  compact  FORTRAN  subroutine. 
This  particular  scheme  has  delivered  good  performance  when  compared  with 
other  methods  (References  4  and  5).  Direct  search  methods  operate  on  the 
basis  of  always  saving  the  most  optimum  point  encountered  as  the  new  "base 
point",  or  point  about  which  further  searches  are  made. 

The  Hooke  and  Jeeves  direct  search  is  unconstrained  in  itself; 
however  as  applied  here  the  problem  is  constrained  through  the  manner 
in  which  the  single-valued  objective  function  is  calculated.  Limits  on  the 
magnitude  of  the  decision  variables,  as  well  as  analytical  relationships 
between  the  decision  variables,  are  imposed  through  the  use  of  individual 
penalties. 

PATSH  operates  by  "moving"  (adjusting)  the  decision  variables 
X1  +  1  »  X1  +  (0.  05)(DEL)  3L1 


where  X1  =  current  decision  variable  set 
Xi+1  =  new  decision  variable  set 
DEL  =  step  size  multiplier 
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These  are  two  results  of  moves.  A  successful  move  produces  a  reduction 
in  the  objective  function  OBJ.  A  move  is  a  failure  when  there  is  no  reduc¬ 
tion  in  OBJ.  Moves  can  be  accomplished  in  one  of  two  ways.  An  explora¬ 
tory  move  consists  of  changing  the  value  of  only  one  decision  variable  and 
evaluating  OBJ.  A  pattern  move  occurs  when  values  of  all  decision  vari¬ 
ables  are  changed  simultaneously  according  to  the  information  derived 
from  exploratory  moves.  During  a  pattern  move,  each  variable  is  changed 
by  an  amount  proportional  to  the  difference  between  its  value  at  the  current 
base  point  and  its  value  at  the  immediately  proceeding  base  point. 

The  logic  flow  of  PATSH  is  presented  in  Figure  8.  PATSH 
begins  the  search  by  calling  the  computational  program  (subroutine  COMP) 
with  the  initial  user- supplied  parameter  set  to  establish  the  initial  base 
point;  this  produces  an  analysis  identical  to  the  first  call  to  COMP  by 
MAIN.  In  the  call  to  the  computational  package,  PATSH  sends  the  current 
parameter  set  to  the  package  and  receives  back  the  objective  function  value 
corresponding  to  that  parameter  set.  After  evaluation  of  the  initial  base 
point,  PATSH  begins  a  series  of  exploratory  moves,  varying  the  value  of 
each  parameter  in  the  following  systematic  manner: 

(1)  Vary  the  parameter  in  the  positive  direction  by  five  percent 
and  evaluate  the  objective  function.  If  the  objective  func¬ 
tion  decreased  in  value  from  the  base  point,  keep  the  param¬ 
eter  change,  save  the  current  total  parameter  set  as  the  new 
base  point,  and  go  to  the  next  parameter. 

(2)  If  the  positive  variation  of  the  parameter  did  not  result  in  a 
reduction  of  the  objective  function,  decrease  the  original 
value  of  the  parameter  by  five  percent  and  evaluate  the 

’  objective  function.  If  the  objective  function  decreases,  a 
new  base  point  is  established;  if  not,  reset  the  parameter 
to  its  original  value  and  go  on  to  the  next  parameter. 

If  the  preceding  exploratory  move  for  this  parameter 
did  not  produce  a  reduction  in  objective  function  when  the 
parameter  was  varied  positively,  but  did  when  it  was 
varied  negatively,  then  the  next  exploratory  move  tries  the 
negative  direction  first  (and  then  the  positive  if  no  improve¬ 
ment  is  seen). 

(3)  When  all  parameters  have  been  varied  one  at  a  time,  eithei* 
a  new  base  point  will  have  been  established,  or  the  original 
base  point  will  be  retained  if  none  of  the  exploratory  moves 
resulted  in  an  improvement.  If  an  improvement  has  been 
achieved,  the  exploratory  moves  have  established  a  pattern  - 
change  the  first  parameter  positive,  do  not  change  the  second 
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parameter,  change  the  third  parameter  negative,  etc,  - 
from  which  a  pattern  move  can  be  taken.  A  pattern  ep 
is  one  in  which  all  parameters  producing  an  improve*  ent 
during  ihe  exploratory  moves  are  varied  simultaneous  . 

If  no  improvement  was  obtained  during  the  exploratory 
moves  (i.  e.  ,  the  previous  base  point  has  been  retained), 
the  step  size  is  reduced  to  one-half  its  current  value  and  the 
exploratory  moves  are  repeated. 

The  pattern  step  may  or  may  not  produce  a  decrease  in  the  objective 
function  over  the  current  base  point.  PATSH  does  not  immediately  reject 
a  pattern  move  that  results  in  an  increase  in  the  objective  function.  Each 
pattern  move  is  followed  by  another  set  of  exploratory  moves,  using  the 
pattern  move  parameter  set  as  the  ’’base"  point.  If  none  of  these  explor¬ 
atory  moves  provides  a  lower  objective  function  value  than  the  base  point 
value  prior  to  the  pattern  move,  the  previous  point  is  retained,  and  a  >et 
of  evoloratory  moves  is  made  about  it.  If  this  set  does  not  produce  a 
r<  ,.;tion  in  OBJ,  the  step  Bize  is  reduced  for  a  new  set  of  exploratory 
moves  about  the  current  base  point.  An  improvement  in  the  objective 
'  »ctior  by  any  means  (exploratory  move  or  pattern  move)  is  always 
-  Gained  as  the  new  base  point.  The  search  is  assumed  to  be  converged 
.en,  through  repeated  efforts  to  obtain  improvements,  the  step  size  is 
uced  from  its  original  value  to  the  minimum  value  specified  by  user 
input  (DELMIN).  Such  a  process  may  appear  to  be  succeeding  by  failing 
to  achieve  any  better  point;  however,  the  final  set  of  exploratory  moves 
clca  iy  demonstrates  no  improvement  in  the  objective  function  by  per¬ 
turbing  all  of  the  parameters  in  either  direction.  This  is  similar  to  evalu- 
at  .  through  a  finite  difference  method,  the  fir9t-order  partial  derivatives 
of  the  objective  function  with  respect  to  the  design  parameters.  Any  error 
in  obtaining  an  optimum  would  be  contained  within  the  minimum  step  size 
used  for  the  final  exploratory  moves. 

When  using  numerical  optimization  techniques,  there  is  always 
concern  over  whether  the  true,  or  global,  optimum  has  been  reached,  >;r 
whether  a  local  optimum  is  the  result.  No  guarantee  exists  that  the  solu¬ 
tion  is  not  a  local  optimum.  The  only  way  to  gain  a  feeling  of  confidence 
in  the  solution  (if  it  is  in  doubt)  is  to  use  different  starting  points  (i.  e. , 
different  initial  (user- supplied)  parameter  sets),  and  to  determine  whether 
or  not  the  same  solution  is  reached  each  time.  The  possibility  of  local 
optima  is  a  function  of  the  problem  to  be  solved.  Some  problems  with 
highly  complex  constraints  may  have  a  number  of  local  optima  while  many 
problems  have  only  one  global  optimum.  Keep  in  mind  that,  even  though  the 
solution  may  be  suspected  to  be  a  local  optimum,  if  all  penalties  are 
zero,  then  the  solution  is  a  valid  design;  some  improvement  in  the  payoff 
parameter  may  be  realized,  and  that  can  be  determined  only  through  start¬ 
ing  the  search  with  a  different  input  set. 
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Performance  Requirement,  Design  Constraint,  and  Operating  Limit 
Satisfaction 


The  optimisation  routine,  PATSH,  operate*  by  minimising  a  single- 
valued  objective  function.  This  single  value  must  reflect  the  pay-off 
quantity  (which  is  multiplied  by  -1  if  maximisation  is  desired)  and  the 
effectiveness  of  the  design  in  meeting  the  performance  requirements, 
design  constraints,  and  operating  limits.  This  has  been  accomplished  by 
incorporating  a  penalty  function  scheme  such  that 

n 

OBJ  *  P  +  Z  T{ 
i=l 


where  OBJ  is  the  single -valued  objective  function  minimised  by  PATSH,  P 
is  the  payoff  quantity,  and  the  Fi  are  individual  penalty  functions  for  each 
of  the  performance  requirements,  design  constraints,  and  operating 
limits  (all  of  which  are  considered  as  constraints  on  the  optimisation  pro¬ 
cess,  and  will  be  referred  to  as  such  for  the  remainder  of  this  discussion). 
Two  basic  types  of  constraints  exist,  inequality  constraints  and  equality 
constraints.  Each  constraint  can  be  expressed  in  one  of  the  following 
forms: 


gj  '  0, 

gt  £ 


For  example,  a  maximum  acceleration  constraint  would  be  expressed  as 


g  =  »  -  a  , 

*a  max  max,  req 

where  am,T  it  the  maximum  acceleration  produced  by  the  design  being 
evaluated  and  am4Jt>  re>;!  is  the  required  maximum.  Clearly,  the  constraint 
is  met  when  ga  $  0,  Similarly,  a  minimum  burnout  velocity  constraint  is 
expressed  as 


g  U  V  -  V 
8V  BO  BO,  req 

This  constraint  is  met  when  g^.  i  0.  Finally,  a  requirement  for  an  exact 
burn  time  would  be  expressed  as 

gtB  =  *B  "  *B,  req 
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and  the  constraint  would  be  satisfied  when 


Penalty  function  values  for  the  three  examples  given  above  would 
He  calculated  as 


Fa  =< 


.  ga  <  0 


ga  sa.  ga>0 


0  »  f?v  *  0 

gy  Sy,  gy  <  0 


Ft 


B 


=< 


gt  =  0 

B 


gt  Stn'  «t  *  0 
lB  B  lB 


The  S&>  Sy,  and  S*  are  scale  factors  used  to  normalize  constraint  viola¬ 
tion  penalties  to  an  appropriate  level  with  respect  to  the  payoff  quantity. 

The  choice  of  this  form  for  the  penalty  functions  provides  a  penalty  value 
that  can  be  scaled  to  relatively  small  values  for  minor  violations  with 
rapidly  increasing  (second  order)  value  for  larger  violations.  Constraint 
enforcement  in  this  manner  can  be  thought  of  as  a  "soft"  constraint  (i.  e.  , 
minor  violations  are  not  totally  excluded  from  the  solution).  Certain  limits 
on  design  parameter  values  are  enforced  as  "hard"  constraints.  An  attempt 
by  the  optimiser  routine  to  specify  a  design  parameter  value  which  violates 
a  "hard"  constraint  results  in  the  specified  value  being  overridden  with  the 
limiting  value  and  the  generation  of  a  penalty  function  proportional  to  the 
attempted  violation.  An  example  of  a  "soft"  constraint  is  the  upper  limit  on 
propellant  web  fraction,  because  a  web  fraction  slightly  greater  than  the 
limit  may  be  acceptable  if  it  produces  greater  improvements  elsewhere. 

An  example  of  a  "hard"  constraint  is  the  length  of  one  part  of  the  motor, 
because  a  length  of  less-than-zero  is  physically  meaningless  (and  can  be 
computationally  misleading). 


Adjustable  Variables 

There  are  36  variables  in  SPOC  which  may  be  adjusted  by  PATSH 
to  obtain  an  optimum  design  (Table  1).  However,  not  all  of  the  decision 
variables  can  be  adjusted  in  any  one  given  problem  because  some  are  pecul¬ 
iar  to  certain  grain  geometries.  The  decision  variables  fall  Into  these 
categories 

o  Propellant  grain  cross-section  dimensions 

o  Propellant  g i  ain  lengths 
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Propellant  ingredient  relative  weights 
Propellant  ballistic  characteristics  (barn  rate 
and  performance  level*  the  latter  as  influenced  by 
ingredient  amounts) 

Noaale  dimensions 

Miscellaneous  (motor  diameter*  case  cylindrical 
wall  thickness) 


PERFORMANCE  REQUIREMENTS 


The  following  performance  parameters  ara  drivan  toward  uaar- 
input  requirement*.  Panaltiaa  ara  calculatad  for  not  mooting  aach  requlre- 
manta.  Dalault  value*  provided  in  tfaa  coda  provant  tha  panaltiaa  from 
being  activated  uhlaaa  the  uaer  chooaaa  to  enforce  the  requirement.  The 
accompanying  parenthetical  expressions  give  tha  appropriate  limit. 

o  Total  impulae  (lower  three- aigma  value  at  low 
temperature) 

o  Total  motor  weight  (maximum  nominal) 

o  Ignition  thruat  (lower  three- aigma  value  at  low  temperature) 

o  Ignition  thruat  (upper  three- aigma  value  at  high  temperature) 

o  Burn  time  (lower-  three -aigma  value  at  high  temperature) 

o  Burn  time  (upper  three- aigma  value  at  low  temperature) 

o  Axial  acceleration  (maximum  nominal  at  high  temperature) 

o  Change  in  velocity  (minimum  nominal  value  at  low  temper¬ 
ature) 

o  Time -to -tar get  (maximum  nominal  value  at  low  temperature) 

o  Impact  (or  termination)  velocity  (minimum  nominal  value  at 
low  temperature) 

Thoae  requirement*  that  are  shown  above  to  apply  to  a  particular 
grain  conditioning  temperature  condition  can  also  be  enforced  with  a  one- 
temperature  problem. 
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DESIGN  CONSTRAINTS  AND 
OPERATING  LIMUS 


Design  constraints  and  operating  limits  that  are  enforced  in  the 
SPOC  arej 

o  Case,  closure  and  noasle  support  thickness  (sufficient 
for  maximum  expected  operating  pressure  plus  safety 
factor) 

c  Case  and  nomale  structure  wall  thickness  (*  manufacturing 
limit) 

o  Noasle  ablative  thickness  (a  that  required  for  char, 
ablation  and  thermal  protection) 

o  Propellant  strain  margin  of  safety  during  low  temperature 
storage  in  both  CP  and  valley  sections  of  grain  (*0) 

o  Propellant  strain  at  low  temperature  ignition  pressurisa- 
tion  (£  input  maximum) 

o  Propellant  web  fraction  (£  maximum  based  on  design 
experience) 

o  Propellant  thickness  under  propellant  valley  (*  manu¬ 
facturing  limit) 

o  Propellant  total  solids  (between  maximum  and  minimum 
limits) 

o  Propellant  burn  rate  and  pressure  exponent  (between 
maximum  and  minimum  limits) 

o  Burn  rate  catalyst  and  fuel  contents  (£  maximum  based 
on  experience) 

o  Combustion  gas  Mach  number  in  port  at  low  temperature 
(nominal  £  maximum  based  on  experience) 

o  Chamber  pressure  at  high  temperature  (nominal  a  maximum 
based  on  experience) 

o  Geometrically  valid  (compatible)  propellant  grain  cross- 
section  dimensions 

o  Lengths  and  thicknesses  greater  than  sero 
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Motor  dimensions  (length  *  maximum,  noasle  exit 
diamatar  a  maximum,  eaaa  alt  opening  radius  •  ncssle 
antranea  radius,  noaala  blast  tuba  length  and  diameter  * 
requirement) 

Geometrically  valid  (compatible)  noasle  dimensions 
Longitudinal  combustion  stability 


PAYOFF  PARAMETERS 


The  PAYOFF  parameter*  from  which  the  user  can  select  one  to  be 
minimized  during  any  given  machine  submission  are 

o  None 

o  Total  motor  cost  (minimize) 

o  Total  motor  weight  (minimize) 

o  Total  impulse  (maximize)^ 

o  Total  impulse- to- total  motor  weight  ratio  (maximize)^  ^ 

o  Burnout  Velocity  (maximize)^ 


(1)  PATSH  will  minimize  the  product  of  minus  one  times  the  value  of  this 
parameter,  which  produces  a  maximization  of  the  parameters. 
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LIMITATIONS  AND  ACCURACY 


The  purpose  of  this  discussion  is  to  summarize  some  general  limi¬ 
tations  of  the  code  and  to  provide  estimates  of  the  accuracy  of  the  results. 

LIMITATIONS 


Some  limitations  on  the  use  of  SPOC  are  inherent  in  the  assumptions 
employed  during  original  development  of  the  analysis  and  simulation 
modules;  these  assumptions  are  given  in  the  discussions  of  the  individual 
modules  and  their  impact  on  a  given  problem  solution  is  best  left  to  the  user. 

Basically,  there  are  no  restrictions  cn  the  size  of  the  motor  which 
may  be  analyzed  with  SPOC.  Small  motors  operating  at  high  pressure  could 
possibly  enter  the  regime  where  thin-wall  pressure  vessel  equations  should 
be  replaced  by  thick- wall  relationships?  it  is  up  to  the  user  to  recognize 
this  situation.  The  cylindrical  section  of  a  motor  employing  elliptical  closures 
(forward  or  aft,  or  both)  cannot  be  reduced  to  zero  length  because  of  how  the 
grain  geometry  is  described  to  the  ballistic  simulation  module;  the  minimum 
length  attainable  is  between  one  and  two  grain  web  thicknesses.  As  for  large 
motors,  there  are  no  restrictions. 

All  volumes  and  concomitant  weights  are  calculated  from  exact  geo¬ 
metric  relationships;  there  are  no  internal  empiricisms  to  estimate  weights. 
Weights  not  amenable  to  direct  calculation  in  an  optimization  code  (e.  g. , 
igniter,  safe-and-arm  device,  wings,  etc.)  are  user-supplied  values. 

Qf  necessity,  some  of  the  analysis  routines  are  somewhat  simplified, 
as  would  b-e  expected  when  operating  in  a  preliminary  design  mode;  how¬ 
ever,  all  analysis  routines  are  industry- accepted  methods. 

(a)  Propellant  strain  is  calculated  under  plane- strain  conditions. 

Thus  end- effects  and  three-dimensional  effects  during  abrupt  configuration 
changes  are  not  acoounted  for. 

(b)  Membrane  stresses  in  the  ellipsoidal  pressure  vessel  closures 
(Type  1)  are  calculated  at  the  motor  centerline  which  provides  a  satisfactory 
estimate  of  the  required  closure  thickness  elsewhere.  Bending  stresses  at 
the  closure-to-cy'iindrical  shell  junction  are  not  considered. 

(c)  Bending  at  the  closure-to-cylindrical  shell  junction  is  considered 
for  the  Type  3  forward  closure  (that  features  a  flat  plate  closure  integral 
with  the  cylindrical  shell)  as  Tong  as  material  response  is  elastic.  Trans¬ 
itions  between  the  cylindrical  shell  and  the  integral  flat  plate  (i.  e. ,  radii  or 
gradually  increasing  cylindrical  wall  thickness  in  the  vicinity  of  the  closure 
are  not  included  in  stress  estimates  or  volume  calculations. 

(d)  The  user  must  input  a  heat-transfer  coefficient  for  each  of  the 
three  nozzle  ablative  materials,  which  means  that  the  coefficient  is  constant 
for  all  flow  conditions  to  which  a  particular  material  is  exposed. 
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There  ere  dimensional  mismatches  at  case-to-elesure  tangent  points 
and  case-to-nozzle  joints  in  order  to  allow  the  user  complete  flexibility 
in  choosing  his  motor  arrangement  and  to  make  the  computations  more 
manageable;  however,  the  results  of  these  mismatches  on  predicted  ballistic 
performance  and  weights  is  thought  to  be  minimal.  Figure  9  shows  the 
potential  mismatch  between  the  liner  inner  surface  at  the  closure-to-caee 
cylindrical  section  interface;  there  are  two  ways  that  this  mismatch  can 
occur,  and  both  are  considered  when  the  grata  outer  dimensions  are  established 
for  the  ballistic  simulation.  Figure  10  shows  the  potential  mismatch  of 
the  pressure  vessel  outer  surface  at  the  closure- to- case  cylindrical  section 
interface.  The  inner  surfaces  of  the  closure  and  cylindrical  section  exactly 
match  at  the  tangent  point i  Then  the  required  closure  thickness  (TCLOF) 
is  calculated  after  the  ballistic  simulation,  and  the  cylindrical  section 
thickness  (TCASE)  is  a  PATSH- adjusted  parameter  that  eventually  is  satis¬ 
factory  for  the  maximum  pressure.  Thus  the  outer  surface  of  the  pressure 
vessel  could  have  a  discontinuity  at  the  tangent  point.  The  thrust  skirt  is 
also  shown  in  Figure  10  to  shew  that  its  mating  surface  is  the  cylindrical 
section  outer  surface.  Obviously,  the  degrees  of  mismatch  shown  in 
Figures  9  and  10  are  greatly  exaggerated  for  clarity;  their  effects  on 
weights  is  negligible. 

Another  mismatch  that  always  occurs  is  shown  in  Figure  11.  The 
case  opening  radius  (RNOZEN)  always  (eventually)  is  equal  to  the  nosale 
entrance  radius,  so  there  is  no  mismatch  there.  However,  the  nozzle 
ablative  and  structural  support  calculations  are  performed  normal  to  the 
internal  surface,  so  that  part  of  the  nosale  coincides  with  the  case  as  shown 
by  the  shaded  area  in  Figure  11;  this  "duplication"  of  volume  provides  an 
allowance  for  the  nozzle  attachment  flange. 

The  trajectory  simulation  employs  a  point-mass  missile  flying  a 
two-dimensional  path  in  the  altitude- range  plane  over  a  flat  earth.  Forces 
modeled  are  reitxicted  to  thrust,  drag  and  weight  (i  .e. ,  lift  is  always  zero), 
and  angle  of  attack  is  always  zero.  The  trajectory  simulation  is  intended  as 
a  supplementary  evaluation  tool  (unless,  of  course,  this  model  accurately 
describes  the  problem  under  consideration). 

a  two-dimensional  plane-strain  model  is  used  to  calculate  propellant 
strain  due  to  low-temperature  storage  and  ignition  pressurization.  Such  a 
model  accurately  describee  the  propellant  behavior  at  a  point  mid-way 
along  the  grain  length  when  the  grain  length-to-diameter  ratio  (L/D)  is 
equal  to  or  greater  than  about  seven.  For  L/D<7,  or  for  locations  near  the 
grain  terminations,  the  plane- strain  models  give  very  conservative  pre¬ 
dictions  because  the  end  effects  (three-dimensional)  that  relieve  the  strain 
are  not  accounted  for  in  SPOC.  Strains  predicted  for  a  propellant  valley  or 
slot  will  also  be  conservative  near  the  ends  or  for  short  slots. 
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Allowance  For  Case 


Figure  11.  Dimensional  Mismatch  at  Case-to- Nozzle  Interface 
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The  propellant  structural  analysis  is  not  conservative  at  the  hinge 
points  of  stress  relief  flaps  and  at  the  transition  between  propellant  slots 
and  CP  regions.  Both  of  these  areas  represent  highly  three-dimensional 
conditions  that  are  uot  amenable  to  preliminary  design  calculations  used  in 
SPOC.  Consequently,  there  is  the  inherent  assumption  that  the  bore  con¬ 
ditions  are  the  critical  locations.  Provisions  have  been  made  to  include 
volume  and  weight  allowances  for  stress  relief  boots  in  ellipsoidal  closures, 
even  though  their  final  configuration  is  dependent  on  more  detailed  analyses. 
The  transition  section  between  slots  and  cylindrical  port  may  require  a 
special  configuration  to  limit  imposed  strains;  another  way  to  achieve  the 
same  results  is  to  specify  about  7  degrees  as  the  angle  on  the  side  of  the 
slot  (ALPHAl)  of  a  finocyl  grain  (Type  3>. 

Thermal  strain  in  the  propellant  due  to  low-temperature  storage  is 
compared  with  design  strain  endurance  (nominal  strain  endurance  reduced 
for  mix-to-mix  variations  and  aging  degradation).  Strain  induced  by  ignition 
pressurization  is  compared  with  a  user- input  maximum  limit.  This  latter 
limit  should  be  derived  from  tests  that  measure  strain  capability  at  rapid 
strain  rate  (to  simulate  ignition  pressurization)  on  test  specimens  condi¬ 
tioned  to  the  design  low  temperature  and  already  strained  to  the  level 
that  will  be  induced  by  low  temperature  storage. 

ACCURACY  OF  CODE 


There  are  three  levels  of  accuracy  tc  consider  in  the  evaluation  of  a 
computer  code.  First,  the  user  must  decide  how  well  the  mathematical 
equations  model  the  reality  of  a  particular  problem.  Second  is  the  compu¬ 
tational  accuracy,  or  how  faithfully  the  programmer  has  carried  out  the 
mathematical  manipulations.  Finally,  and  totally  under  the  control  of  the 
user,  is  the  accuracy  of  the  input  data,  ftily  the  first  two  levels  will  be 
discussed  here. 

Accuracy  of  the  mathematical  models  is  paramount  in  the  overall 
accuracy  of  a  code.  The  several  analysis  and  simulation  modules  are  dis¬ 
cussed  separately  in  the  following  list: 
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Module 


Estimated  Accuracy  of  Model 


Ballistic  simulation 


Weight  estimates 


Propellant  theoretical 
characteristics 

Combustion  stability 


Combustion  efficiency 


Motor  costs 


T  r  a  j  ector  imulation 


Propellant  structural 
analysis 


Press  '  ■«  vessel  «•  > 
tural  analysis 


13%  total  impulse 
19%  maximum  pressure 

General  qualitative  assessment  based  on  experience 
12% 

General  qualitative  assessment  based  on  experience 

Essentially  error  free.  Uses  NASA- Lewis 
thermochemical  analysis  (Reference  11). 

Based  on  AFRPL  Standard  Stability  Code 
(Reference  19). 

Based  on  AFRPL  Solid  Propellant  Prediction 
Code  (Reference  12). 

Based  on  Tri-Services  Rocket  Motor  Trade-off 
Study  for  steel  cases  (Reference  18). 

Estimated  to  be  very  high,  provided  the  problem 
is  adequately  described  by  the  model.  See 
discussion  above. 

Strain  calculation  "very  accurate"  in  center  of 
motor  with  LD  >7  (probably  within  10%).  For 
location  near  ends  of  long  motor  or  for  L/D  <  7, 
calculated  strains  are  conservative,  with  degree 
of  conservatism  depending  on  problem. 

Estimated  to  be  conservative  by  approximately 
15%. 


The  computational  accuracy  of  the  code  is  extremely  high.  Iteration 
schemes  in  the  bal)*;<tics  simulation  and  grain  subroutines  require  conver¬ 
gence  to  within  0.0  >r  less.  The  trajectory  simulation  uses  an  industry- 
accepted  technique.  Thus  it  is  felt  that  the  mathematical  models  have  jeen 
faithfully  computed. 
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MOTOR  GEOMETRY 


This  section  of  Volume  I  provides  more  detailed  discussion  of  the 
motor  geometry  aspects  of  SPOC,  as  differentiated  from  descriptions  of  the 
various  analysis  routines.  As  with  much  of  the  manual,  some  information 
given  here  applies  to  other  parts  of  the  code,  and  so  there  is  some  dupli¬ 
cation. 


The  user  m?1'  choose  one  from  each  of  three  forward  closure  arrange¬ 
ments,  two  aft  closure  arrangements  and  six  nozzle  configurations  to 
describe  the  basic  motor  geometry.  All  of  these,  in  any  combination,  may 
be  used  with  any  of  the  five  propellant  grain  configurations,  except  that  a 
conocyl  grain  can  be  employed  only  with  an  ellipsoidal  (Type  1)  forward 
closure  (simply  because  of  the  definition  of  a  conocyl  grain). 

FORWARD  CLOSURE 


There  are  two  basic  forward  closure  arrangements:  an  ellipsoid  and 
a  flat  plate.  The  ellipsoidal  closure  is  Type  1.  The  flat  plate  closure  pro¬ 
vides  two  types)  one  has  the  closure  secured  in  the  case  by  a  separate 
retaining  ring  (Type  2)  and  the  other  has  the  closure  integral  with  the  cylin¬ 
drical  portion  of  the  case  (Type  3). 

Type  1  Forward  Closure 

The  controlling  interface  between  the  Type  1  forward  closure  and 
the  cylindrical  portion  of  the  case  is  at  the  tangent  point  of  the  ellipsoid  and 
the  cylinder  on  the  inner  surface  of  the  pressure  vessel  (Figure  12);  all 
else  derives  from  this.  The  semi-major  axis  is  found  from 

B2F  =  RMOTOR  -  TCASE  (1) 

Then  using  the  input  ellipse  ratio,  the  semi-minor  axis  of  the  controlling 
interface  surface  is 


A2F  =  B2F/BE1 A2F  (2) 

These  two  parameters  then  are  used  to  define  the  .temi-major  and  semi-minor 
axes  of  other  ellipsoids  using  input  and  calculated  material  thickness  (see 
Figure  12), 
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FORWARD  CLOSURE  TYPE  1 


ire  Type  1  (Ellipsoid) 


A1F  »  A2F  +  TCLOF 

(3) 

A3F  =  A2F  TINFMX 

(4) 

A4F  *  ASF  -  TSRBF 

(5) 

A5F  *  A4F  -  TLNR 

(6) 

B1F  =  B2F  +  TCLOF 

(7) 

B3F  »  B2F  -  TINFMN 

(8) 

B4F = B3F 

B5F  =  B4F  -  TLNR 

(10) 

It  is  assumed  that  the  subsequent  bodies  thus  described  are  ellipsoids,  which 
is  not  so  mathematically,  but  in  reality  the  difference  is  very  small  for 
any  situation  where  the  thicknesses  used  in  Eqs.  (3)  •  (10)  are  small  compared 
with  the  semi-minor  axis. 

The  thicknesses  TINFMX,  TSRBF,  and  TLNR  are  defined  at  the 
igniter  opening  RIGN;  however,  to  establish  the  ellipse  parameters,  these 
thicknesses  are  assumed  to  act  at  the  motor  centerline.  This  approach  does 
not  introduce  significant  error  because  the  igniter  opening  is  usually  small 
with  respect  to  the  motor  sices. 

The  closure  thickness,  TCLOF,  is  calculated  after  the  ballistic  simu¬ 
lation  is  performed  so  that  the  stress  analysis  may  use  the  design  pressure 
for  the  current  pass  through  COMP.  The  other  ellipsoidal  surfaces  are 
defined  prior  to  the  ballistic  simulation  because  they  are  needed  to  describe 
the  propellant  outer  surface  in  the  closure  for  the  ballistic  simulation  routine. 
Radii  (RECH)  and  concomitant  lengths  (HECH)  measured  forward  from  the 
first  plane  that  describes  the  propellant  grain  are  internally  calculated 
(Figure  13).  This  first  plane  is  positioned  aft  of  the  closure/case  tangent 
point  a  distance  of  TAUMXF,  which  is  the  maximum  distance  burned  in  any 
given  propellant  grain. 

It  can  be  noted  on  Figure  13  that  there  is  no  attempt  to  achieve  a 
perfect  blend  at  the  tangent  point  between  identical  materials  located  in  the 
closure  and  in  the  cylindrical  portion  of  the  case.  The  liner  inner  surface 
does  not  necessarily  match  at  the  tangent  point,  etc.  The  error  in  propellant 
and  insulation  volumes  should  be  small  for  any  reasonable  combination  of 
input  thicknesses,  and  the  computation  routine  is  greatly  simplified. 

Because  of  some  rigidly  enforced  geometric  rules  in  the  ballistic  simula¬ 
tion  routine,  this  mismatch  of  the  liner  inner  surface  makes  a  difference 
where  the  first  coordinate  (RECH(l),  HECH(l))  is  placed  (Figure  14). 

The  innermost  surface  of  the  propellant  (e.  g. ,  the  box*e  of  a  CP 
grain)  must  intersect  the  liner  inner  surface  on  a  Type  l  closure.  There  are 
no  provisions  for  a  forward  opening  larger  than  the  propellant  inner  surface. 
Thus  the  control  input  to  inhibit  burning  surface  in  the  forward  closure 
(FWDINH)  has  no  meaning  for  a  Type  1  closure,  and  so  it  is  bypassed. 


47 


HECH(4)rec^ 


Type  2  and  Type  3  Forward  Clotuf 


Type  2  and  Type  3  closure*  (Figure  IS)  are  identical  with  regard  to 
the  propellant- closure  arrangement;  they  are  differentiated  only  in  how  the 
closure  i*  secured  in  the  case  to  form  the  pressure  vessel.  The  controlling 
interface  between  the  Type  2  and  Type  3  closures  and  the  cylindrical  portion 
of  the  case  is  the  outer  edge  of  the  flat  plate,  the  inner  surface  of  the  case, 
and  the  aft-facing  surface  of  the  flat  plate. 

The  forward  face  of  the  propellant  grain  is  always  perpendicular  to 
the  motor  centerline;  it  can  be  inhibited  (FWDINH  ■  T)  or  not  (FWDINH  =  F). 
A  space  can  be  provided  between  the  closure  insulation  and  propellant  with 
the  input  LGAPF.  The  inputs  to  simulate  a  grain  cast  against  and  bonded 
to  the  closure  are  LGAPF  =  0  and  FWDINH  =  T. 

Although  not  shown  in  Figure  15,  the  liner  and  insulation  under  the 
grain  extends  forward  across  LGAPF  until  it  contacts  the  closure  insulation. 

The  first  plane  describing  the  grain  for  the  ballistic  simulation  is 
positioned  a  distance  TAUMXF  aft  of  the  forward,  propellant  surface.  Only 
two  radii  and  lengths  are  needed  to  describe  the  propellant  forward  of 
Plane  1.  HECH(l)  and  HECH(2)  are  equal  to  TAUMXF.  RECH(l)  is  equal 
to  the  propellant  outside  radius  at  Plane  1.  If  the  forward  propellant  face 
is  inhibited,  RECH(2)  is  zero;  if  it  is  not  inhibited,  RECH(2)  is  equal  to 
RECH(l). 

AFT  CLOSURE 


There  are  two  basic  aft  closure  arrangements.  Type  1  is  an  ellipsoid 
with  an  opening  for  nozzle  attachment;  Type  2  causes  the  closure  to  be  formed 
by  the  entrance  section  of  the  nozzle. 

Type  1  Aft  Closure 

The  controlling  interface  between  the  Type  1  aft  closure  and  the 
cylindrical  portion  of  the  case  is  at  the  tangent  point  of  the  ellipsoid  and  the 
cylinder  on  the  inner  surface  of  the  pressure  vessel  (Figure  16.  Definition 
of  ellipse  parameters  proceeds  as  with  the  ellipsoid  forward  closure  with 
one  significant  exception.  The  thicknesses  TINAMX,  TSRBA,  and  TLNR 
are  measured  normal  to  the  local  surface  at  the  nozzle-to-case  interface, 
RNOZEN;  these  are  then  used  to  establish  the  ellipses. 

The  closure  thickness.  TCLOA,  is  calculated  after  the  ballistic 
simulation  is  performed.  The  other  ellipsoid  surfaces  are  defined  prior  to 
the  ballistic  simulation.  Radii  (RECN)  and  concomitant  lengths  (HECN)  are 
measured  aft  from  the  last  plane  that  describes  the  propellant  grain.  This 
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HE  AD  END  TYPES  2  AND  3 


Figure  15.  Forward  Closure  Type  2  and  Type  3  (Flat  Plate) 


AFT  CLOSURE  TYPE  1 


Aft  Closure  Type  I  (Ellip: 


last  plane  is  positioned  forward  of  the  aftmost  propellant  surface  a  distance 
of  TAUMXA  (Figure  17). 

Note  on  Figure  17  that  the  axial  length  of  the  aft  closure,  LCLOA, 
extends  aft  from  the  case** to- closure  tangent  point  to  the  interface  of  the 
stress  relief  boot  and  the  closure  insulation.  This  point  also  serves  as  the 
interface  between  the  case  and  nozzle  (at  the  radius  RNOZEN). 

As  for  the  Type  1  forward  closure,  the  mismatch  between  liner  inner 
surfaces  of  the  case  and  closure  is  considered  through  the  internal  genera^- 
tion  of  appropriate  radii  and  lengths  (Figures  17  and  18). 

Aft-facing  burning  surface  perpendicular  to  the  motor  centerline  is 
formed  when  the  aft  case  opening,  RNOZEN,  is  larger  than  the  innermost 
surface  of  the  propellant  (e.  g. ,  the  bore  of  a  CP  grain).  This  perpendicular 
surface  can  be  inhibited  (AFTINH  =  T)  or  not  (AFTINH  =  F). 

Type  2  Aft  Closure 


An  aft  closure  Type  2  is,  in  actuality,  not  a  closure  per  se  (Figure 
19);  instead  the  pressure  vessel  is  dosed  by  the  entrance  section  of  the 
nozzle.  There  is  no  specific  controlling  interface  between  the  Type  2  aft 
closure  and  the  case.  The  controlling  interface  with  the  nozzle  is  at  the 
outside  propellant  radius  at  Plane  14  (RFA14).  In  other  words,  the  case 
"opening**1  for  the  nozzle,  RNOZEN,  is  set  equal  to  RFA14.  Length  of  the 
aft  closure  is  zero,  which  means  that  the  nozzle  is  positioned  at  the  pro¬ 
pellant  aft  face. 

The  aft  face  of  the  grain  is  always  perpendicular  to  the  motor  center- 
line;  it  can  be  inhibited  (AFTINH  =  T)_or  not  (AFTINH  =  F).  The  laBt  plane 
describing  the  propellant  grain  for  the  ballistic  simulation  is  positioned  a 
distance  TAUMXA  forward  of  the  aft  propellant  surface.  Only  two  radii  and 
lengths  are  needed  to  describe  the  propellant  aft  of  Plane  14.  HECN(l) 
and  HECN(2)  are  equal  to  TAUMXA.  RECN(l)  is  equal  to  the  propellant  out¬ 
side  radius  at  Plane  14.  If  the  aft  propellant  face  is  inhibited,  RECN(2)  is 
zero;  if  there  is  no  inhibiting,  RECN(2)  is  equal  to  RECN{1). 

Thickness  of  the  aft  skirt  ic  TCASE  for  a  Type  2  aft  closure  (for 
purposes  of  skirt  weight  calculations).  It  was  reasoned  that  for  most  appli¬ 
cations  with  this  aft  end  arrangement  the  skirt  would  be  an  extension  of  the 
case. 

NOZZLE 


Six  nozzle  configurations  were  identified  for  user  s 
20  through  25).  A  review  of  recent  designs  (References  6 


election  (Figures 
and  7)  and 
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Figure  20.  Dimensional  Inputs  for  Nozzle  Type 
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the  CPIA  Manual  (Reference  8)  indicates  these  configurations  should  describe 
most  applications  for  which  SPOC  is  intendsd. 

Nossle  Type  1  and  Nossle  Type  2  (Figures  20  and  21)  are  basically 
the  same,  except  Type  2  provides  for  a  one-piece  throat  insert  in  a  conical 
seat*  Type  1  is  more  akin  to  all-plabtic  noszles  where  low  weight  is  impor¬ 
tant,  although  the  computation  routine  allows  any  material  to  be  used  for 
structural  support.  Nossle  Type  3  (Figure  22)  is  derived  from  smaller 
tactical  motors  where  the  support  structure  is  an  extension  of  the  case. 
However,  it  is  allowable  to  attach  the  Type  3  nossle  to  Type  1  aft  closure. 
Nossle  Type  4  (Figure  23)  was  established  to  provide  missile  equipment 
volume  around  the  outside  of  the  nossle  and  to  have  its  throat  located  at 
the  forward  end  of  the  nossle  section.  An  identical  external  envelope  is 
found  with  Nossle  Type  5,  but  the  nossle  throat  is  located  at  the  aft  end  of 
the  blast  tube.  Nozsle  Type  6  still  provides  equipment  volume  at  the  aft  end 
of  the  motor,  but  now  the  envelope  allows  an  exit  cone  to  be  attached  to  the 
arrangement  of  the  Type  5  nossle. 

All  nossle  types  have  the  option  of  either  conical  or  contoured 
expansion  sections.  In  the  former,  inputs  of  ALFA  =  ALFAEX  describes 
the  conical  section  to  the  nozzle  routine,  and  a  flag  is  set  by  the  user 
(CONTUR  =  F)  as  part  of  the  motor  definition  inputs.  The  user  specifies  a 
contoured  expansion  section  by  CONTUR  =  T  and  ALFAEX  <  ALFA.  When 
a  contoured  section  is  specified,  the  code  internally  chooses  either  an 
elliptical,  hyperbolic  or  parabolic  profile  on  the  basis  of  minimum  length. 
Neither  of  these  contours  will  be  identical  to  one  determined  by  precise  gas 
dynamic  analyses,  but  the  estimates  of  weights  and  lengths  are  sufficiently 
accurate  for  the  purposes  of  SPOC. 

As  described  earlier,  the  nozzle  entrance  radius,  RNl,  is  eventually 
made  equal  to  the  aft  case  opening,  RNOZEN;  at  the  beginning  of  an  optimi¬ 
zation  problem  there  may  be  points  where  RNl  and  RNOZEN  are  not  equal, 
but  m  appropriate  penalty  is  calculated  to  force  them  together. 

LENGTHS 

Motor  length  is  summed  as  shown  in  Figures  26  and  27,  These 
two  illustrations  show  a  Type  5  (CP)  grain,  but  the  technique  is  consistent 
for  all  grains.  Individual  lengths  between  Plane  1  and  Plane  14  are  unique 
for  a  particular  grain  configuration.  Length  components  forward  of  Plane  1 
and  aft  of  Plane  14  are  identical  for  all  grain  configurations,  as  is  the  defin¬ 
ition  of  LCASE. 

Figure  28  illustrates  the  details  of  the  case-nozzle  interface  for  a 
Type  1  aft  closure  and  how  resultant  lengths  are  defined.  The  nozzle  entrance 
is  always  placed  at  the  junction  of  the  stress  relief  boot  and  closure  internal 
insulation.  Thus  LCLOA  is  measured  from  that  point  forward  to  the  case- 
closure  tangent  point,  and  LNOZ  is  measured  from  that  point  aft  to  the  nozzle 
exit  plane. 
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Figure  27.  Typical  Length  Summation  for  Forwa 

Aft  Closure  Type  2 


Figure  28.  Case-Nozzle  Interface  for  Type  1  Aft  Closure 


66 


ANALYSIS  ROUTINES 


This  section  of  Volume  I  describes  the  analyses  performed  within 
SPOC  to  determine  how  the  design  currently  being  proposed  by  PATSH 
satisfies  performance  requirements,  design  constraints  and  operating  limits. 
Propellant  grain  configuration  dimensions  are  checked  for  geometric  validity 
and  are  converted  into  the  language  of  the  ballistic  simulation  subprogram. 
Propellant  ingredient  weight  fractions  are  checked  for  compatibility  with 
interacting  limits,  normalized  to  100%  total  weight  fraction  and  then  sent  to 
the  thermochemical  subprogram.  There  ballistic  characteristics  and 
parameters  for  impulse  efficiency  and  combustion  stability  predictions 
are  calculated.  Impulse  efficiency  is  predicted  with  the  AFRPL  Solid 
Propulsion  Prediction  (SPP)  code  (Reference  12)  and  combustion  stability 
with  the  AFRPL  Standard  Stability  Prediction  (SSP)  code  (Reference  19). 

A  ballistic  simulation  is  performed,  including  the  effects  of  erosive  burning 
and  mass  addition  and  using  a  geometrically  rigorous  two-dimensional  grain 
regression  technique.  Once  the  ballistic  simulation  is  completed,  results 
from  it  are  used  for  pressure  vessel  and  propellant  structural  analyses, 
a  trajectory  simulation,  cost  estimate,  and  weight  calculations,  and  com¬ 
bustion  stability  prediction. 

Some  analyses  are  not  performed  ut.less  specified  by  the  user: 
thermochemical,  impulse  efficiency,  propellant  structural,  trajectory, 
cost,  and  combustion  stability.  Even  when  the  thermochemical  analyses  is 
called,  it  is  bypassed  unlesB  the  pressure  or  the  nozzle  expansion  has 
changed  at  least  5%  between  the  two  previous  analyses,  or  if  propellant 
formulation  has  changed. 

Equations  describing  the  various  analyses  are  numbered  within 
each  report  section.  All  references  are  listed  at  the  end  of  this  volume. 
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BALLISTIC  SIMULATION 


The  ballistic  simulation  routine  requires  inputs  of  the  initial  propel¬ 
lant  configuration,  propellant  ballistic  and  burning  rate  properties,  nozzle 
geometry  and  various  control  parameters.  The  primary  computed  outputs 
are  chamber  pressure,  thrust,  mass  flow  rate,  pressure-time  integral, 
impulse,  and  propellant  weight  at  selected  values  of  burning  time. 

The  initial  propellant  geometry  is  described  at  a  series  of  planes 
(called  "direct  input"  planes)  positioned  along  the  grain  in  accordance  with 
a  pre-programmed  scheme.  Fourteen  planes  are  used  for  all  grain  config¬ 
urations;  they  are  located  as  shown  in  the  illustrations  describing  each  con¬ 
figuration.  Planes  not  shown  on  the  illustrations  are  called  "interpolated 
input"  planes  because  they  are  positioned  equidistance  between  and  their 
dimensions  are  internally  derived  from  the  two  adjacent  "direct  input"  planes 

The  gas  dynamic  solution  is  based  on  the  assumption  of  equilibrium 
(i.  e. ,  rate  of  mass  stored  is  negligible)  at  each  time  point.  Other  assump¬ 
tions  are: 

a.  Combustion  products  behave  as  perfect  gases. 

b.  Flow  processes  are  isentropic, 

c.  Propellant  burning  surface  regresses  only  in  a  radial 
direction,  except  on  the  forward»facing  and  afUfaclng  sur¬ 
faces  at  the  forward  and  aft  grain  terminations,  respectively 

d.  Burning  rate  is  calculated  at  each  plane  and  varies 
linearly  between  planes. 

e.  Ends  of  propellant  grain  do  not  experience  erosive 
burning. 

f.  AU  interpolations  are  linear. 

g.  Chamber  is  already  filled  and  ail  of  the  initial  burning 
surface  is  ignited  at  the  initial  (time  =  zero)  calculation. 

h.  Contribution  of  igniter  mass  flow  to  initial  pressure 
is  ignored. 

i.  Two- phase  mixture  effects  and  other  losses  are  accounted 
for  through  an  empirical  impulse  efficiency  factor. 

Operation  is  initiated  with  a  solution  for  incremental  burning  surface 
areas  and  propellant  volumes  throughout  the  motor.  Then  a  trial  value  of 
stagnation  chamber  pressure  at  the  head- end  of  the  motor  is  assumed. 

Based  on  the  assumed  pressure,  incremental  values  of  burning  rate,  gas 
flow  rate,  temperature,  velocity,  specific  weight,  and  Mach  number  are 


calculated  at  the  adjacent  station  in  the  downstream  direction.  These  calcu¬ 
lations  include  the  effects  of  both  mass  addition  pressure  drop  and  erosive 
burning  effects  on  propellant  burning  rate.  This  solution  is  repeated  for 
each  successive  station  along  the  grain  until  a  solution  is  obtained  for  gas 
properties  at  the  aft  end  of  the  grain.  Gan  discharge  rate  through  the  noz¬ 
zle  is  computed  from  the  nozzle  end  stagnation  pressure,  using  the  appro¬ 
priate  solution  for  either  sonic  or  subsonic  flow.  The  gas  discharge  rate 
is  then  compared  to  the  gas  generation  rate.  If  these  two  values  do  not 
agree  within  a  prescribed  tolerance,  the  original  trial  value  of  head- end 
pressure  is  adjusted  to  a  new  value  and  the  entire  solution  is  repeated 
until  the  required  agreement  is  obtained.  Motor  thrust  is  then  computed, 
and  all  inputs  are  printed.  Burning  time  is  incrementally  advanced  and  the 
thickness  burned  at  each  plane  is  calculated,  based  on  the  respective  burn¬ 
ing  rate  values  previously  calculated.  New  incremental  values  of  burning 
surface  and  volume  are  computed,  and  the  entire  ballistic  analysis  is 
repeated.  This  process  is  continued  until  all  of  the  propellant  has  been 
consumed,  or  until  a  specified  time  or  pressure  level  is  obtained. 


Ballistic  simulations  can  be  performed  at  either  one  or  two  grain 
soak  temperatures.  If  simulations  are  run  at  two  temperatures,  results 
at  the  high  temperature  are  used  to  calculate 


a.  Upper  three- sigma  ignition  thrust, 

b.  Lower  three- sigma  burn  time, 

c.  Maximum  head-end  pressure, 

d.  Pressure  vessel  design  pressures  (MEOP,  yield, 
ultimate). 


and  results  at  the  low  temperature  are  used  to  calculate 


a.  Lower  three- sigma  total  impulse 

b.  Lower  three- sigma  ignition  thrust 

c.  Upper  three- sigma  burn  time 

d.  Maximum  port  Mach  number 

e.  Upper  three-sigma  ignition  thrust  for  propellant  structural 

analysis 

f.  Burn  time  for  nozzle  thermal  analysis 


If  a  simulation  is  run  at  only  one  temperature,  all  of  the  above  are  calculated 
at  that  single  temperature. 


Each  grain  configuration  has  its  own  individual  subroutine  (SETUPl, 
SETUP? ,  etc. )  that  generates  the  initial  propellant  geometry  description 
from  the  fewest  possible  inputs.  The  different  configurations  are  described 
in  detail  in  the  following  eections.  All  grain  subroutines  use  another  sub¬ 
routine  (CLOS)  that  performs  the  time  function  for  the  different  pressure 
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vessel  closure  configurations  that  may  be  selected.  In  addition,  CLOS  cal¬ 
culates  the  weights  of  inert  components  associated  with  the  closures,  as 
does  each  of  the  "SETUP"  routines  for  their  respective  grain  configurations. 
Propellant  surface  regression  and  internal  ballistics  are  calculated  at  four 
stations  in  the  clcsures,  as  shown  on  Figures  29  and  30. 

Burning  rate  is  calculated  with  the  Vielle  relationship 

YW 

RATE  =  (BRSF)  (A)  P  (1) 

The  user  supplies  RB70  (the  rate  at  1000  psia,  70*F)  and  XN.  The 
code  uses  conventional  temperature  sensitivity  coefficients  to  adjust  the 
coefficient  "A"  to  temperature  extremes. 

Erosive  burning  can  be  considered  through  selection  of  one  of  the 
combinations  shown  on  Table  2  (References  9  and  10). 

When  either  the  second  or  third  options  of  Table  2  are  selected,  the 
burn  rate  becomes 

RATE  =  (BRSFMAjP3™  (M/MCRIT)XM  (2) 

or, 

MPFVP 

RATE  (MPCOEFJ(MP)  (3) 

where  MP  is  the  product  of  the  local  Mach  number  and  static  pressure. 

The  relation  producing  the  largest  burn  rate  at  a  given  location  and  time 
is  used  to  calculate  the  internal  gas  dynamics  (either  Eq  (1),  (2),  or  (3)). 

Nozzle  throat  ablation  is  modeled  by 

RE  =  (KREl)(PON)  (4) 

where  PON  is  the  nozzle-end  stagnation  pressure.  The  nozzle  thermal 
analysis  subroutine  also  calculates  an  ablation  profile  for  the  entire  nozzl* 
internal  surface,  but  that  is  not  used  in  the  ballistic  simulation. 

Impulse  efficiency  is  a  user  input,  or  it  may  be  estimated  by  the 
subroutine  IMPEFF  which  uses  the  AFRPL  SPP  "empirical"  model  (Reference 
12).  Certain  inputs  to  the  SPP  impulse  efficiency  calculation  must  come 
from  a  thermochemistry  analysis  of  a  particular  propellant  formulation; 
if  a  formulation  is  uot  input  to  the  code,  the  user  must  supply  values  for 
the  specified  parameters. 
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Station  4 


Forward  Closure  Type  1 


Figure  29.  Location  of  Stations  for  Ballistic  Simulation  in  Head- end  of 
Grain 
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TABLE  2 


EROSIVE  BURNING  RATE  COMBINATIONS 


AUTOEB 

MCRIT 

XM 

MPCOEF 

MPEXP 

Description 

0.0 

0.0 

0.  0 

0.0 

0.0 

No  erosive  burning 

0.0 

x<*> 

X 

X 

X 

Erosive  burning  with 
user  values  for  MCRIT, 
XM,  MPCOEF  and 
MPEXP 

1.0 

0.0 

0.  0 

X 

X 

Erosive  burning  with 

internally  calculated 
value*  for  MCRIT  and 
XM  (»ee  Note  (b)  )  and 
user  values  for 
MPCOEF  and  MPEXP 
(see  Note  (c)  ) 


(a)  "X"  indicates  user- input  or  default  value. 

(b)  Saderholm  model.  Reference  9  and  10.  MCRIT  ~  Mach  number  corre¬ 
sponding  to  250  fps.  XM  calculated  from  the  relationship. 


XM  =  In  0. 
curve  fit  to  S 


/p°.  74 

06763  ~ - 
l  RAT  El 


0.  4948 


which  is  a  curve  fit  to  Saderholm  data  (XM  *  0,  0  and  RATE  =  burn  rate 
without  erosive  burning,  P  =  local  static  pressure,  when  local  Mach 
number  >  MCRIT. 

(c)  Default  values  MPCOEF  *  0.  0093  and  MPEXP  =  0.  71. 
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Nozzle  divergence  loa a  is  calculated  according  to 

f  ALFA  +  ALFAEXl 

LAMDA  =  1  -  COS  2  j  (5) 

where 

ALFA  =  Half-angle  (deg)  at  entrance  to  nozzle  expansion 

section 

ALFAEX  =  Half-angle  (deg)  at  exit  of  nozzle  expansion 

section 

This  term  LAMDA  is  included  in  the  calculation  of  thrust  coefficient  within 
the  ballistic  simulation  module,  unless  the  SPP  Impulse  efficiency  option  is 
selected.  Because  SPP  efficiency  includes  divergence  losses,  LAMDA 
per  Eq  (5)  is  internally  set  equal  to  1.0. 

Certain  other  gas  dynamic  parameters  may  be  input  by  the  user  or 
calculated  internally  if  a  propellant  formulation  is  supplied. 

The  ballistic  simulation  uses  inputs  from  the  following  namelists: 

a.  BALLST 

b.  GRAIN  x  (x=l,...  5) 

c.  NOZGEC 

d.  IN G AMT  (if  thermochemical  analysis  performed) 

e.  Card  that  names  ingredients 


PROPELLANT  GRAIN  CONFIGURATIONS 


Basic  describing  dimensions  for  the  propellant  grains  are  furnished 
by  the  user,  who  can  then  specify  that  certain  cf  them  be  adjusted  by  PATSR 
during  the  optimisation  process.  Regression  of  the  burning  surface  is 
mathematically  exact,  and  that  means  that  the  dimensions  furnished  to  the 
ballistic  simulation  module  must  obey  strict  rules.  If  the  rules  are  not 
obeyed,  the  burning  surface  regression  cannot  proceed  and  the  run  is  aborted. 
Abnormal  terminations  are  unwanted  in  a  pattern  search  optimization  process. 

Therefore,  a  subroutine  was  formulated  for  each  of  the  different  grain 
configurations  to: 

(a)  check  the  geometric  validity  of  the  incoming  dimension  set, 

(b)  adjust  c  'rtain  dimensions  to  obtain  geometric  validity  (if 
required),  and  calculate  associated  penalties 

(c)  derive  from  the  incoming  dimensions  those  other  dimensions 
required  by  the  ballistic  simulation  module 

(d)  compare  the  dimensions  with  design  constraints 

(e)  calculate  volume  and  weights  of  inert  components  associated 
with  the  cylindrical  portion  of  the  pressure  vessel 

(f)  calculate  grain  dimensions  that  are  needed  in  subsequent 
analyses 

The  grain  geometries  are  described  in  the; ballistic  simulation 
routine  by  a  series  of  planes  oriented  perpendicular  to  the  motor  centerline 
and  positioned  along  the  length  of  the  grain  at  appropriate  locations.  Four¬ 
teen  planes  are  used  to  describe  all  grains.  Some  are  located  where  there 
are  changes  in  either  the  internal  or  external  configuration  of  the  propellant: 
these  are  knowh  as  "direct  input"  planes,  and  it  is  their  dimension  sets 
that  are  treated  in  the  aforementioned  geometric  validation  analysis.  The 
regaining  planes  are  located  equidistance  between  adjacent  "direct  input" 
planes:  they  are  known  as  "interpolated  inputs"  because  their  dimensions 
are  derived  from  linear  interpolations  from  the  adjacent  direct  input  planes. 
The  planes  that  are  indicated  on  the  illustrations  accompanying  this  dis¬ 
cussion  are  the  direct  input  planes  established  for  each  of  the  grain  con¬ 
figurations. 

Each  grain  type  has  its  own  particular  arrangement  of  internal 
insulation.  Details  of  the  insulation  are  given  in  Volume  II  of  this  Manual 
and  will  not  be  repeated  here. 

Dimensions  that  can  be  adjusted  during  the  optimization  process 
are  discussed  in  detail  in  subsequent  paragraphs.  The  capabilities  for 
adjustment  are  not  the  same  for  all  the  grains. 
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(a)  Lengths  of  the  various  grain  segments  can  be  adjusted  for  all 
grains  during  the  jptlmization  process. 

(b)  For  the  Type  1  (Star)  grain,  all  major  dimensions  describing 
the  cross-section  are  either  held  constant  at  the  initial  input, 
or  they  are  all  varied  during  the  optimization  process. 

(c)  For  the  Type  2  (Wagon  Wheel)  grain,  all  dimensions  describing 
the  cross-section  are  either  held  constant  at  the  initial  input,  or 
they  are  all  varied  during  the  optimization  process. 

(d)  For  the  Type  3  (finocyl)  grain,  the  major  dimensions  describing 
the  cross-section  can  be  individually  held  constant  or  allowed 
to  vary  during  the  optimization  process. 

(e)  For  the  Type  4  (conocyl)  grain,  the  major  dimensions  describing 
the  cross-section  can  be  individually  held  constant  or  allowed 

to  vary  during  the  optimization  process. 

(f)  For  the  Type  5  (CP)  grain,  the  dimensions  describing  the 
cross-section  can  be  individually  held  constant  or  allowed 
to  vary  during  the  optimization  process. 

Basic  arrangements  of  the  closures  are  ellipsoidal  and  flat  plate. 

Either  of  these  arrangements  can  be  used  with  any  of  the  five  grain  configur¬ 
ations  except  that  the  Type  4  (conocyl)  grain  can  be  used  only  with  the 
Type  1  (ellipsoidal)  forward  closure;  this  is  because  the  definition  of  a  conocyl 
grain  positions  the  slot  adjacent  to  the  forward  closure  and  the  geometric 
grain  regression  routine  was  established  for  this  arrangement. 

Otherwise,  the  treatment  of  the  closures  is  the  same  for  all  grains. 
The  illustrations  giving  the  plane  locations  show  that  any  variation  in  grain 
dimensions  along  the  length  of  the  motor  (e.  g.  port  radius,  star  height) 
does  not  continue  forward  of  Plane  1  or  aft  of  Plane  14.  In  other  words,  the 
configuration  that  exists  at  Plane  14  is  "projected"  into  the  aft  closure,  even 
though  in  the  actual  motor  the  variation  in  port  radius  (for  example)  would 
probably  continue  until  the  end  of  the  grain  in  the  aft  closure.  This  technique 
slightly  overpredicts  the  volume  of  propellant.  However,  Plane  1  and 
Plane  14  are  positioned  only  a  short  distance  of  TAUMXF  and  TAUMXA, 
respectively,  from  the  grain  ends  where  TAUMXF  and  TAUMXA  are  tne 
maximum  distances  burned  at  Plane  1  and  Plane  14,  respectively. 

Type  1  (Star)  Grain 


Type  1  grain  is  a  standard  star  configuration  (Figure  31).  Loca¬ 
tions  of  direct  input  planes  are  shown  in  Figures  32  and  33  for  the  two 
closure  types.  Dimensions  shown  in  blocks  on  Figure  31  are  for  Plane  1, 
which  are  held  constant  to  Plane  7;  then  the  height  of  the  star  tip  tapers 
from  LSA1  at  Plane  7  to  LSA14  at  Plane  14.  All  other  dimensions  are 
constant  from  Plane  7  to  Plane  14. 


Notes 

(1)  Line  AB  is  perpendicular  to  Line  AC. 

(2)  Dimensions  in  blocks  arc  input;  others  are  output 


Figure  31. 


Type  1  (Star)  Grain  Configuration 


Figure  33.  Location  of  Direct  Input  Planes  for  Type  1  (Star)  Grain  with  Type 


Subroutine  SETUP  1  performs  the  validation  checks  described  above 
for  the  Type  1  grain.  The  logic  flow  for  the  cross-section  dimension 
checking  is  shown  in  Figure  34. 

There  are  two  categories  of  defining  variables.  First  are  those 
that  are  fixed  (insofar  as  SETUP1  is  concerned)  and  are  never  adjusted.  If 
these  data  are  found  to  be  unacceptable,  the  run  is  terminated  with  an 
appropriate  message. 

RFAl 


R2A1 

THETA 

W 

WFLIM 


Propellant  outside  radius.  Is  actually 
adjusted  through  changes  in  motor  diameter 
and  case  thickness,  but  SETUP1  has  no 
control  over  these. 

User-selected  input 
180/NSLOTS.  See  Figure  31. 

Minimum  clearance  between  adjacent 
star  tips 

User-selected  input  of  maximum  allowable 
web  fraction 


Second  are  those  dimensions  that  are  varied  by  PATSH  during  the  optimi¬ 
zation  process.  If  these  data  are  found  to  be  unacceptable,  they  are  adjusted 
in  accordance  with  a  hierarchy  described  below. 


TAUW1 

LSA1 

R5A1 

ALPHA1 


Web  thickness 

Star  point  height  above  web 

Fillet  radius  between  star  point  and  web 

Included  half-angle  of  star  point 


The  first  group  of  tests  involve  the  fixed  inputs. 


RFAl  >  0  (1) 

R2A1  *0  (2) 

THETA  >0  i3) 

W  20 

(R2A1  +  W/2)  <  (RFAl)Sin{ THETA)  (4) 

0  <  WFLIM  <1.0  (5) 


Failure  of  any  of  the  above  tests  results  in  a  run  termination. 
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The  second  group  of  testa  ares 

0  <  TAUWl  <  (WFLIM)(RFA1)  (6) 

R2A1  <  LSAl  <  LSMAX  (?) 

where  LSMAX  =  RFA1  -  TAUW1  +  R2A1  -  (R2A1  +  W /2) /Sin(THET  A)  (8) 

Failure  of  Eq.  (6)  causes  TAUW1  to  be  set  to  the  limit.  Failure  of  Eq.  (7) 
causes  LSAl  to  be  set  to  the  limit.  Penalties  are  calculated  in  both 
instances. 

At  this  point,  only  ALPHAl  and  R5A1  of  the  four  adjustable  variables 
remain  to  be  tested.  The  third  group  of  tests  is  concerned  with  ALPHAl. 
The  decision  logic  is  shown  on  Figure  35.  The  variable  LSSTAR  is  the 
maximum  star  tip  height  possible  when  ALPHAl  =  THETA,  for  the  situation 
where  RFAl  Cos  (THETA)  <  (RFAl  -  TAUW1)  and  RFAl,  TAUWl,  R2A1 
and  THETA  have  been  previously  accepted  (Figure  36).  The  following 
possibilities  can  exist: 

If  LSAl  >  LSSTAR 

If  LSAl  =  LSSTAR 

If  LSAl  <  LSSTAR 


Maximum  ALPHAl  must  be  less 
than  THETA 

Maximum  ALPHAl  equals  THETA 
and  R5A1  is  unique  at  input  value 
Maximum  ALPHAl  equals  THETA 
and  multiple  solutions  exist  for 
RFAl 


When  LSAl  is  greater  than  LSSTAR  (the  only  situation  where  maximum 
ALPHAl  i  THETA),  maximum  ALPHAl  is  determined  from  the  criteria 
shown  on  Figure  37.  Note  from  Figure  35  that  minimum  ALPHAl  is 
always  zero.  If  ALPHAl  in  the  current  dimension  set  is  less  than  zero  or 
greater  than  maximum  ALPHAl,  it  is  set  to  the  limit  and  an  appropriate 
penalty  calculated. 

The  fourth  group  of  tests  is  concerned  with  R5A1,  the  last  dimension 
needed  to  define  the  cross-section  geometry.  A  portion  of  the  R5A1  logic 
diagram  is  given  in  Figure  38.  The  first  step  is  to  solve  for  the  angle 
AFSTAR  (Figure  39),  where  R5A1  can  be  zero,  subject  to  the  criteria 
shown  on  Figure  39.  Thus, 

ALPHAl  <  AFSTAR  Minimum  R5A1  =  0 

ALPHAl  >  AFSTAR  Minimum  R5A1  >  0 

Therefore,  the  next  step  is  to  determine  minimum  R5A1  when 
ALPHAl  >  AFSTAR,  and  it  is  shown  in  Figure  40.  The  final  step  of  this 
group  of  tests  is  to  determine  maximum  R5A1.  The  solution  for  maximum 
R5A1  depends  on  the  regime  in  which  ALPHAl  is  located;  this  regime  is 
defined  by  the  angle  DELTA,  which  is  shown  on  Figure  41.  DELTA 
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IF  RFAl  cob  (THETAi<  \RFA1  -  TAUW1)  ?'F  RFAl  cos  (THETA)  a 

I  (RFAl  -  TAUW1) 


SOLVE  FOR  LSSTAR  (SEE  FIGURE  006) 


IF  LSAl  >  LSSTAR  IF  LSA1  *  LSSTAR 

I _ 

SOLVE  FOR  MAX  ALPHAl  MAX  ALPHA  1  =  THETA 

(SEE  FIGURE  007)  MIN  ALPHAl  =  0.0 

I 

MIN  ALPHAl  =  0.  0 


FORCE  ALPHAl  TO  LIE  WITHIN  RANGE  SHOWN 
MIN  ALPHA  *  ALPHAl  «  MAX  ALPHAl 


Figure  35.  Testing  of  ALPHAl  for  Type  1  Grain 
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Figure  36. 


Definition  of  LSSTAR,  Type  1  Grain 


Limiting  Criteria: 


A  line  from  the  center  of  R5A1  (which  is  Point  A) 
and  perpendicular  to  line  BD  must  not  intersect 
line  CE  at  a  point  further  from  Point  C  than 
Point  E. 

Subject  to:  RFAl,  TAUWl,  THETA,  R2A1,  LSAl  already  accepted 


Figure  37.  Definition  of  Maximum  ALPHAl  when  LSAl  >LSSTAR,  Type  1  Grain 
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SOLVE  FOR  AFSTAR 
(SEE  FIGURE  39) 


IF  ALPHA  1  *  AFSTAR 
MIN  R5A1  =  0 


IF  ALPHA1  >  AFSTAR 


SOLVE  FOR  MIN  R5AI 
(SEE  FIGURE  40) 


J 


I - 1 - 

IF  R5A1  <  MIN  R5AI  IF  R5Al  =  MIN  R5A1 


RESET  R5A1 
TO  EQUAL  MIN  R5A1 


DATA  HAS  PASSED 
INPUT  CHECKS 


IF  R5A1  >  MIN  R5A1 


GO  TO  SOLUTION  FOR 
R5.MAX 


I 


Figure  38. 


Logic  Diagram  tor  AFSTAR,  Type  1  Grain 
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Limiting  Criteria:  A  line  from  Point  B  and  {Perpendicular  to  line  AB 

cannot  inter eect  line  CAD  at  a  point  further  from 
Point  C  than  Point  D 

Subject  to:  RFAl,  TAUW1,  THETA,  R2A1,  LSAl  already  accepted 


Figure  39. 


Definition  of  AFSTAR,  Type  1  Grain 


Limiting  Criteria:  A  line  from  the  center  of  R5A1  (which  ie  Point  A) 

and  perpendicular  to  line  BE  muit  not  intersect 
line  CD  at  a  point  further  from  Point  C  than 
Point  D. 


Figure  40. 


Minimum  RFAl  when  ALPHAl  >  AFSTAR,  Type  1  Grain 


(illustrated  in  part  B  of  Figure  41)  is  the  unique  value  of  ALPHAl  at  which 
the  center  of  R5A1  radius  lies  on  the  line  DE  and  a  line  tangent  to  the 
R2A1  and  R5A1  arcs  is  perpendicular  to  a  line  connecting  their  centers. 

For  ALPHAl  >  DELTA,  (part  A  of  Figure  41),  the  center  of  R5A1  radius 
lies  oh  the  line  DE,  but  a  line  connecting  the  centers  of  R2A1  and  R5A1  is 
not  perpendicular  to  a  lire  tangent  to  the  R2A1  and  R5A1  arcs.  For 
ALPHAS  DELTA  (part  C  of  Figure  41),  the  center  of  R5A4  radius  lies 
inside  line  DE,  and  a  line  tangent  to  the  R2A1  and  R5A1  arcs  iB  perpendicular 
to  a  line  connecting  their  centers.  A  test  is  made  to  determine  in  which 
regime  ALPHAl  is  located,  and  appropriate  geometric  relationships  are 
solved  to  find  maximum  R5A1.  Then,  the  final  test  is  to  insure  that 
R5MIN  <  R5A1  <  R5MAX.  If  R5A1  is  outside  the  limits,  it  is  set  to  the  limits 
and  a  penalty  is  calculated. 

Now  that  all  dimensions  for  Plane  1  are  acceptable,  star  point  height 
at  Plane  14  (LSA14)  is  tested  to  insure  R2A1  *  LSA14  sLSAl,  where  R2Al 
is  the  same  minimum  height  requirement  at  Plane  1. 

The  final  tasks  of  subroutine  SETUP  1  are  to  insure  that  grain 
lengths  are  greater  than  zero,  calculate  inert  weights,  calculate  data  for 
use  in  propellant  structural  analysis,  translate  the  grain  dimensions  into 
language  required  by  the  ballistic  simulation  module,  and  other  miscellaneous 
tasks.  The  dimension  L3  (Figure  42)  is  derived  for  use  in  propellant  struc¬ 
tural  analysis,  through  the  following  expression; 

LITD  =  (L3+R5)(2)  (9) 

so  that  LITD  becomes  the  distance  across  the  propellant  valley.  The  distance 
T  max  shown  on  Figure  42  is  the  maximum  distance  burned  for  the  star  grain 
ancl  ie  used  to  establish  the  positions  of  the  first  and  last  planes  (that  describe 
the  grain  to  the  ballistic  simulation  module)  with  respect  to  the  forward  and 
aft  closures  (e.  g.  ,  see  Figure  26). 


91 


|— LIT  D/2-] 

Dimensions  for  Propellant  Structural  Analy 
Type  1  Grain 
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Type  2  (Wagon  Wheel)  Grain 


Type  2  grain  is  a  double-web  wagon-wheel  configuration  (Figure 
43).  Locatior*  of  direct  input  planes  are  shown  in  Figures  44  and  45  for  the 
two  closure  types.  Dimensions  shown  in  blocks  on  Figure  43  are  for  Plane 
1,  which  are  held  constant  to  Plane  7;  then,  the  height  of  the  star  point  can 
taper  from  LSAl  at  Plane  7  to  LSA14  at  Plane  14.  All  other  dimensions  are 
constant  from  Plane  7  to  Plane  14  as  long  as  they  are  geometrically  compatible 
with  the  reduced  LSA14.  Basic  dimensions  describing  the  grain  cross-section 
(TAUWl,  LSAl,  R5A1)  can  all  be  varied  during  the  optimization  process, 
or  they  all  may  be  held  constant  at  the  initial  input. 

Subroutine  SETUP2  performs  the  validation  checks  described  above 
for  the  Type  2  grain.  The  logic  flow  for  the  cross-section  dimension  check¬ 
ing  is  shown  in  Figure  46  for  all  variaoles  being  adjusted  and  in  Figure  47 
for  all  variables  held  constant.  Although  two  different  flow  charts  are  shown, 
the  two  routines  have  been  merged  in  the  coda. 

There  are  four  categories  of  defining  variables.  First  are  those  that 
are  fixed  (insofar  as  SFTUP2  is  concerned)  and  are  never  adjusted.  If  these 
data  are  found  unacceptable,  the  run  is  terminated  v/ith  an  appropriate  mes¬ 
sage. 


RFA1  >  0  Propellant  outside  radius.  Is  actually 

adjusted  through  motor  diameter  and 
case  thickness,  but  3ETUP2  has  no 
control  over  these. 

THETA  >  0  180/NSLOTS 

W 02  MIN  0  Half  the  minimum  clearance  between 
adjacent  star  points 

Second  are  those  input  data  that  are  accepted  if  possible,  but  SETUP2 
will  adjust  them  if  required  to  avoid  a  run  termination, 

WFLIM  User- selected  input  of  maximum 

allowable  web  fraction. 

This  input  is  accepted  if  possible.  If  not,  it  is  adjusted  to  be  compatible 
with  the  remaining  dimensions  in  order  to  prevent  run  termination.  It  is 
reset  to  the  appropriate  values  with  no  penalties  being  calculated.  A 
geometrically  maximum  possible  wab  fraction  (independent  of  the  user  input) 
exists  when  TAUWl  increases  to  the  point  where  R5A1  is  driven  to  zero  and 
W02  is  equal  to  W02MIN  (Figure  48);  this  condition  results  when  the  length 
LC  (shown  in  Figure  43)  is  zero.  At  this  condition,  web  thickness  TAUWl 
is  tl.e  largest  possible  and 


^TAUMAX  shown  in  Figure  48  is  not  the  maximum  distance  burned  used  to 
position  ballistic  planes;  the  latter  has  the  same  definition  as  in  Figure  42. 
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gure  44.  Location  of  Direct  Input  Planes,  Type  2  Grain 


PUu  Number 


HECN(1 ) 


Figure  46.  (Continued) 


-W  02  MIN 


Figure  48.  Maximum  Web  Thickness  TAUMAX  For  Wagon  Wheel 
(Type  2)  Grain 


YTT 


(9) 


WFMAX  - 


TAUMAX 

RFA1 


WFLIM  is  reset  to  WFMAX,  or  held  at  its  input  value,  whichever  is  smaller. 


Third  are  the  variables  adjusted  by  PATSH  during  the  optimisation 
process.  These  data  are  accepted  as  input,  if  possible;  but,  if  they  are  not 
geometrically  acceptable,  they  are  reset  to  appropriate  limits  and  penalties 
are  calculated. 


TAUW1  Weu  thickness  at  Plane  1 

LSA1  Star-point  height  above  web  at  Plane  1 

R5A1  Fillet  radius  between  star  point  and 

web  at  Plane  1 

R2A1  Corner  radius  at  top  of  star  point 

Tho  fourth  and  final  category  contains  those  dimensions  that  are 
calculate  j  in  SETUP2  for  transmittal  to  the  ballistic  simulation  module. 

They  are  based  on  all  other  dimensions. 

LA  Length  of  star  point  (See  Figure  43) 

AL1A1  Included  half- angle  or>  star  point  (See 

Figure  43)  at  Plane  1 

The  first  step  in  the  dimension  check  is  to  test  RFAl,  THETA,  and 
W02MIN  as  shown  above.  If  these  tests  are  passed,  W02  is  set  equal  to 
WOZMIN  and  TAUMAX  (as  defined  in  Figure  48)  is  found.  Knowing  TAUMAX, 
the  maximum  web  fraction  WFMAX  is  found  from  Eq  (9)  and  compared  to  the 
input  WFLIM.  WFLIM  is  set  equal  to  WFMAX  if  the  latter  is  smaller  than 
the  input  value.  Next,  web  thickness  TAUWI  is  forced  to  be  within  the  range 

0  s  TAUWI  i  (WFLIM)(RFA1) 

If  PATSH  has  adjusted  TAUWI  outside  these  limits,  TAUWI  is  reset  to  the 
limit  and  a  penalty  calculated. 

The  maximum  possible  corner  radius  R2A1  is  TAUWI  in  order  to 
maintain  the  star  point  with  the  same  web  thickness  as  TAUWI.  So  the  next 
test  is 


0  s  R2A1  £  TAUWI 

If  PATSH  has  adjusted  R2A1  outside  these  limits,  R2A1  is  reset  to  the  limit 
and  a  penalty  calculated. 
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The  fourth  atep  ia  to  aolve  for  the  length  LA  (shown  in  Figure  43). 


TAUW1  -  R2A1 
ten  (THETA) 


(10) 


Note  that  when  R2A1  =  TAUWl,  LA  =  0, 


Now  that  R2A1  end  TAUWl  ere  acceptable  (or  have  been  react  to 
acceptable  valuea),  W02  ia  made  equal  to  W02MIN  ao  that  the  maximum  pos- 
aible  atar  point  height  (LSMAX)  car*  be  calculated  (aee  Figure  49).  The 
minimum  poaaible  atar  point  height  (LSMIN)  ia  alao  calculated  (aee  Figure 
50);  LSMIN  occura  for  a  given  TAUWl  when  R5A1  and  LC  are  aero.  Then 
the  incoming  atar  point  height  (LSA1)  ia  forced  to  be  within  the  range 

LSMIN  <  LSAl  &  LSMAX 


If  PATSH  haa  adjusted  LSAl  outaide  theae  limits,  LSAl  ia  reaet  to  the  limit 
and  a  penalty  ia  calculated. 

The  final  teat  ia  to  aaaure  that  the  fillet  radius  R5A1  is  between 
maximum  and  minimum  limits  that  are  determined  from  the  particular 
dimensional  combination.  Minimum  R5A1  is  zero;  although  small  radii 
are  critical,  strains  calculated  in  the  propellant  structural  analysis  will 
serve  bo  limit  the  minimum  fillet  radius.  Maximum  R5A1  can  be  established 
by  one  of  two  situations.  The  center  of  the  radius  R5  must  lie  within  the  pie- 
shaped  segment  of  symmetry  or  on  its  boundary.  In  the  first  situation  (Figure 
51),  R5MAX  occura  when  the  length  LC  equals  zero  and  the  center  of  R5  ia 
within  the  segment  of  symmetry.  In  the  aecond  situation  (Figure  52),  the 
center  of  R5  ia  on  the  boundary  of  the  segment  of  symmetry  and  LC  i  0. 
Whichever  situation  prevails,  R5A1  ia  forced  to  lie  within  the  range 

0  s  R5A1  s  R5MAX 


If  PATSH  has  adjusted  R5A1  outside  theae  limits,  R5A1  is  reset  to  the  limit 
and  a  penalty  is  calculated. 


102 


Now  that  all  dimensions  for  Plane  1  are  acceptable,  star- point  height 
at  Plane  14  (I..SA14)  is  tested  to  insure  LSMIN  *  LSA14  £  LSA1.  If  LSA14  is 
anything  other  than  LSAl,  the  checks  described  above  for  LSA.1  are  per¬ 
formed  again  for  LSA14,  and  R5A14  (different  from  R5A1)  is  calculated,  if 
necessary. 

The  final  tasks  of  subroutine  SETUP2  are  to  insure  that  grain  lengths 
are  greater  than  zero,  calculate  inert  weights,  calculate  data  for  use  in  pro¬ 
pellant  structural  analysis,  translate  the  grain  dimensions  into  language 
required  by  the  ballistic  simulation  module,  and  perform  other  miscellaneous 
tasks. 


Type  3  (Finocyl)  Grain 

Type  3  grain  is  a  finocyl  configuration  (Figure  53)  with  the  longi¬ 
tudinal  slots  located  in  the  forward  end  of  the  grain.  Locations  of  direct 
input  planes  are  shown  in  Figures  54  and  55  for  the  two  closure  types 
Dimensions  shown  in  blocks  on  Figure  53  are  for  Plane  1,  which  are  held 
constant  to  Plane  4.  Between  Plane  4  and  Plane  7,  the  slot  depth  radius  R5 
decreases  until  at  Plane  7  it  equals  the  bore  radius  R2.  Slot  fillet  radii  at 
Planes  5  and  6  (R4A5  and  R4A6)  are  calculated  internally  to  be  tangent  to  the 
slot  sides  formed  by  angle  ALPHAl.  The  circular  port  radius  (R21)  is  held 
constant  between  Planes  7/8  and  10/11,  and  then  it  can  expand  at  Plane  14  to 
control  gas  velocity. 

The  dimensions  varied  by  PATSH  during  the  optimization  process 
are  R5A1,  R4A1,  R2A1,  and  R2A14.  Subroutine  SETUP3  performs  the 
validation  checks.  The  logic  flow  for  the  cross-section  dimension  checking 
is  shown  in  Figure  56. 

A  restriction  on  the  minimum  size  of  the  bore  radius  R2A1  is  a  ueer- 
input  of  RIGN,  that  provides  a  clearance  for  the  igniter  and  gas  flow  passage 
around  it.  The  slot  fillet  radius  R4A1  ha*  a  lower  limit  R4MIN  also  supplied 
by  the  user.  An  indirect  restriction  on  the  minimum  R4A1  is  the  strain 
level  calculated  in  the  propellant  structural  analysis  (which  is  then  compared 
to  an  allowable  strain).  The  user  also  supplies  a  MINWEB  dimension  .which 
sets  the  maximum  R4A1  for  a  given  RFA1.  The  minimum  R5A1  is  R2A1  plus 
R4A1.  The  slot  side  angle  ALPHAl  is  constant  at  the  user  input  value,  but 
it  is  compared  with  ALPHMX  for  every  incoming  dimension  set;  if  ALPHAl 
>  ALPHMX,  it  is  set  equal  to  ALPHMX  for  the  current  dimension  set.  It  is 
then  reset  back  to  input  ALPHAl  for  the  next  dimension  set. 

Type  4  (Conocyl)  Grain 

Type  4  grain  is  a  conocyl  configuration  with  the  transverse  slot 
located  in  a  Type  1  (ellipsoidal)  forward  closure  (Figure  57).  Locations  of 
the  direct  input  planes  are  given  in  Figures  58  and  59  for  the  two  closure 
cypes.  Dimensions  shown  in  blocks  on  Figure  57  are  input  and  result  in  a 
description  of  the  conocyl  grain. 
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Figure  56 


Continued  (Part  6  of  6) 


Figure  57.  Nomenclature  for  Head  End  of  Grain  Type  4 


Direct  Input  Plane 


A  geometrically  exact  burn  surface  regression  model  for  the  forward 
portion  of  the  grain  produces  a  history  of  burning  surface  and  port  area  (the 
latter  at  Plane  2/3)  versus  distance  burned.  This  burning  surface  is  con¬ 
verted  to  a  perimeter  versus  distance  burned  (needed  by  ballistic  simulation 
module)  by  dividing  the  surface  area  by  the  length  L2T1,  where  L2T1  = 
TAUMXF  +  A5F,  and  A5F  is  the  semi-minor  axis  of  the  liner  inside  surface 
(see  Figure  12).  The  bore  radius  R2  is  held  constant  to  Plane  9,  and  then 
it  can  expand  at  Plane  14  to  control  gas  velocity. 

Subroutine  SETUP4  performs  the  validation  checks  for  the  Type  4 
grain.  The  logic  flow  for  the  forward  segment  checking  is  shown  in  Figure 
60.  These  checks  are  illustrated  on  Figure  6l.  At  the  forward  end,  the 
slot  tip  is  tested  to  assure  it  is  located  within  the  bounds  shown  in  Figure  61, 
Part  K,  and  that  the  slot  angle  ZED  is  45°<ZED<90°  (Part  E).  Other 
checks  «ire  performed  at  the  aft  end  so  that  R2A9<R2A14<RFA14  (Part  F 
and  Part  H)  and  that  the  nozzle  entrance  radius  RNOZEN  is  properly  esta¬ 
blished  (Part  G).  The  other  checks  assure  RIGN  is  dimensionally  compatible 
with  the  remainder  of  the  problem  and  that  port  radius  R2A3  is  witihin  limits. 

The  final  tasks  of  subroutine  SETUP4  are  to  insure  that  grain  lengths 
are  greater  than  zero,  calculate  inert  weights,  translate  the  grain  dimensions 
into  language  required  by  the  ballistic  simulation  module,  and  perform  other 
miscellaneous  takks. 

Type  5  (CP)  Grain 

Type  5  grain  has  a  cylindrical  port  (CP)  for  its  entire  length.  Direct 
input  planes  ar«  shown  in  Figures  62  and  63  for  the  two  closure  types. 
Cross-sectional  dimensions  adjusted  in  the  optimization  process  are  the  port, 
radii  at  Planes  1,  5  and  14  (R2A1,  R2A5,  R2A14,  respectively).  Tne  dimen¬ 
sional  checks  are  relatively  simple. 

R2A1  *  R2A5 
R2A14  *  R2A5 

Other  cheekt  determine  that  R2MIN^  R2  <  RF,  where  R2MIN  is  established 
by  a  web  fraction  limit. 

Subroutine  SETUPS  performs  these  verifications.  As  with  the  other 
grains,  SETUP5  also  insures  that  lengths  aed  greater  than  zero,  calculates 
inert  weights  associated  with  the  cylindrical  portion  of  the  pressure  vessel, 
translates  the  grain  dimensions  into  language  required  by  the  ballistic  simu- 
latior  'dule,  and  performs  other  miscellaneous  tasks. 
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Figure  62.  Direct  Input  Planes  for  Type  5  Grain,  Closure 


propellant 


Propellant  ballistic  and  gas  dynamic  properties  can  be  supplied  to 
SPOC  by  several  means. 


Option  1  allows  evaluation  of  a  single  propellant  whose  characteristics 
are  known.  It  would  be  employed  when  the  user  wants  to  design  a  motor  that 
contains  a  specific  propellant.  The  flag  to  set  for  Option  1  is  PROPIN  =  T 
and  two  other  flags  (described  beiow)  are  allowed  to  default  to  FORMAD  =  F 
and  FORMIN  =  F.  The  user  must  furnish  propellant  ballistic  characteristics: 


CSTR70 

DELP 

GAMAC 

RGAS 

RB70 

XN 


Characteristic  velocity  (nozzle  end)  at  70 °F 
(ft/sec) 

Propellant  cured  density  (Ibm/cu  in) 

Ratio  of  specific  heats  in  chamber 
Gas  donstant  of  combustion  products  in  chamber 
(ft-lbf/lbm-  °R) 

Propellant  burn  rate  at  70 *F,  1.000  paia  (in/ sec) 
Pressure  exponent  in  the  rate  model 
RATE  =  A*P**  N 


If  the  problem  requires  ballistic  simulation  at  two  grain  temperatures,  addi¬ 
tional  user- supplied  inputs  are; 

PIK  Temperature  coefficient  of  pressure  (per  °F) 

MC  Temperature  coefficient  of  characteristic 

velocity  (per  #F) 

If  the  SPP  impulse  efficiency  model  (Reference  12)  is  specified  by  the  user 
(SPPETA  =  T),  the  user  must  supply  additional  information  for  this  propel¬ 
lant  option;: 


IVAC 

IVACF 

MOLCND 


Vacuum  specific  impulse,  shifting  equilibrium 
at  motor  pressure  and  expansion  ratio 
(Ibm-  sec/lbm) 

Vacuum  specific  impulse,  frozen  equilibrium 
at  motor  pressure  and  expansion  ratio 
(lbf- sec/lbm) 

Mole  fraction  of  condensed  speci&r  (moles  per 
100  gms  of  mixture) 


If  a  combustion  stability  analysis  is  to  be  performed  (FSTAB  =  T),  the  user 
must  supply  more  information  for  this  option: 


SONVEL  Sonic  velocity  in  chamber,  stagnation  conditions 

(ft/sec) 
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Option  2  will  permit  the  uier  to  supply  the  formulation  of  a  given  pro¬ 
pellant  and  SPOC  will  calculate  much  of  the  information  listed  above.  This 
option  is  also  employed  when  the  user  wants  to  design  a  motor  that  contains 
a  specific  propellant,  as  did  Option  1.  The  flag  to  set  for  Option  2  is 
FORMIN  =  T,  and  two  other  flags  are  allowed  to  default  to  FORMAD  =  F 
(described  below)  and  PROPIN  =  F.  In  this  option,  the  thermochemistry 
module  (TCHEM),  is  called  only  on  the  first  pass  through  COMP  and  these 
results  are  used  throughout  the  optimisation  process.  Chamber  pressure 
and  nozzle  expansion  ratio  are  needed  in  TCHEM;  on  the  first  pass  through 
COMP  the  chamber  pressure  that  is  used  is  an  input  value,  PC  (defaults  to 
1000  psia),  and  nozzle  expansion  ratio  is  derived  from  initial  nozzle  dimen¬ 
sions.  Data  that  must  still  be  supplied  by  the  user  are  RB70,  XN,  PIK,  and 
MC.  All  other  data  are  produced  by  the  thermochemistry  module,  provided 
another  condition  is  satisfied:  IVAC,  IVACF  and  MOLCND  are  not  calculated 
unless  SPPETA  =  T. 

Option  3  is  employed  when  the  propellant  ingredient  weight  fractions 
are  adjusted  as  part  of  the  optimization  process.  The  flag  to  set  for  Option  3 
is  FORMAD  =  T,  and  two  other  flags  are  allowed  to  default  to  FORMIN  =  F 
and  PROPIN  =  F.  In  this  option,  the  thermochemistry  module  (TCHEM)  is 
entered  on  every  pass  through  COMP  (subject  to  the  conditions  discussed 
below)  to  provide  the  data  listed  above.  Chamber  pressure  and  nozzle  expan¬ 
sion  ratio  are  needed  in  TCHEM.  On  the  first  pass  through  COMP,  the 
chamber  pressure  that  is  used  is  an  input  value,  PC,  and  nozzle  expansion 
ratio  i»  derived  from  the  initial  nozzle  dimensions.  On  subsequent  passes 
through  COMP,  the  chamber  pressure  is  the  average  pressure  from  the 
preceeding  simulation,  and  the  expansion  ratio  is  still  based  on  the  current 
nozzle  dimension  set;  if  the  problem  being  solved  calls  for  ballistic  simulation 
at  two  temperatures,  the  average  pressure  used  is  from  the  high-temperature 
simulation,  As  the  optimum  design  is  approached,  the  difference  between  the 
pressures  from  the  current  simulation  and  the  just  completed  simulation  will 
disappear.  The  thermochemistry  module  is  not  called  for  a  new  analysis 
unless  the  chamber  pressure  or  the  nozzle  expansion  ratio  on  the  last  pass 
through  COMP  is  at  least  5%  different  from  that  on  the  next-to-last  pass  or 
the  ingredient  weight  fractions  have  been  adjusted. 

TCHEM  is  the  executive  subroutine  for  the  entire  thermochemistry 
module.  It  call?  subroutines  LIQUID  and  NORMAL  (described  below)  and 
then  MAINCO,  which  is  the  executive  subroutine  for  the  thermochemistry 
analysis  itself  (once  a  valid  formulation  has  been  provided  to  it). 

Theoretical  density  calculated  by  the  thermochemical  analysis  routine 
is  based  on  ingredient  density  and  relative  amounts.  This  value  is  multiplied 
by  0.  985  to  obtain  a  "cured"  density  to  account  for  polymerization  and  cool¬ 
down  from  cure  temperature. 

Propellant  Formulation 

The  propellant  formulation  is  adjusted  as  part  of  the  optimization 
process  through  changes  in  weight  fractions  of  related  ingredients.  Two 
steps  are  required  in  this  process.  First,  the  current  formulation  must  be 
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checked  against  user-supplied  and  physical  limits  and  then  normalized  bo  that 
the  total  weight  fraction  is  1.0.  Second,  this  verified  data  set  is  sent  to  the 
thermochemistry  subroutine  to  have  'the  various  properties  calculated. 

The  first  task  performed  by  subroutine  TCHEM  1b  (through  subroutine 
LIQUID)  to  verify  that  the  incoming  formulation  is  physically  possible  and  that 
it  adheres  to  all  limits  (Figure  64).  In  thiB  effort,  it  checks  the  compatibi¬ 
lity  of  liquid  constituents  and  total  solids  content.  For  example,  if  all  liquid 
ingredients  have  their  weight  fractions  fixed  by  the  uaer  at  input  values,  that 
sets  the  total  solids  level;  if  that  exceeds  a  user- supplied  maximum  total 
solids  content,  a  basic  incompatibility  exist;  and  subroutine  LIQUID  so 
informs  the  user. 

The  second  task  performed  by  subroutine  TCHEM  is  (through  sub¬ 
routine  NORMAL)  to  normalize  the  incoming  propellant  ingredient  weight 
fractions  so  that  the  total  is  1.0  (Figure  65).  One  or  more  weight  fractions 
can  be  changed  by  PATSH  during  the  optimization  process.  Subroutine 
NORMAL  takes  these  quantities,  (knowing  the  ingredients  to  be  adjusted  by 
PATSH  and  the  limits  calculated  by  subroutine  LIQUID)  to  specify  a  formu¬ 
lation  whose  weight  fractions  total  1.  0  and  that  are  within  user-  or  internal¬ 
generated  limits.  It  calculates  penalties  when  incoming  weight  fractions  must 
be  changed  to  conform  to  limits.  Only  those  ingredients  being  adjusted  bv 
PATSH  are  changed  in  subroutine  NORMAL  to  normalize  the  formulation. 

The  normalization  process  is  an  iterative  one  (Figure  65). 

Within  the  thermochemistry  analysis  module  itself  is  another  normali¬ 
zation  routine  that  will  adjust  all  ingredient  weight  fractions  to  total  1.0. 
regardless  of  t  :  intent  of  the  user.  This  situation  should  not  be  encountered 
unless  the  formulation  is  input  with  the  weight  fraction  total  not  equal  to  1.  0 
and  flags  set  to  hold  all  ingredients  at  their  input  weight  fractions 
(FORMIN  =  T).  Then,  the  second  normalization  process  will  change  all 
amounts  so  that  the  total  is  equal  to  1.0. 

The  next  task  performed  by  subroutine  TCHEM  is  to  set  up  the  pro¬ 
pellant  ingredient  data  prior  to  calling  MAINCO  (the  entry  to  the  ther;  10- 
chemical  analysis).  Up  to  this  point,  oxidizer  has  been  considered  as  up  to 
two  independent  materials  with  up  to  three  independent  particle  sizes  of 
each;  now,  these  are  combined  into  total  amounts  of  each  cf  the  two  oxidizers. 

Thermochemistry 

Propellant  ballistic  and  gas  dynamic  characteristics  are  calculated 
internally  in  SPOC  (when  so  specified)  using  a  version  of  the  NASA-Lewis 
code  TRAN72  (Reference  01).  In  its  original  form,  this  code  calculates 
thermodynamic  and  transport  properties  of  complex  mixtures,  chemical 
equilibrium  for  assigned  thermodynamic  states,  and  theoretical  rocket  per¬ 
formance  for  both  frozen  and  equilibrium  compositions  during  expansion  for 
condensed  and  gaseous  species.  Ingredient  and  specie  thermodynamic  and 
transport  properties  are  obtained  from  JANNAF  tables  that  are  periodically 
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Figure  65. 


Logic  Flow  Chart  for  Subroutine  NORMAL 


updated.  Execution  time  was  reduced  by  streamlining  the  code  so  that  it 
performs  only  those  calculations  needed  by  SPOC,  produces  no  printout, 
and  considers  only  those  combustion  products  that  will  be  produced  by  a 
limited  group  of  propellant  constituents.  This  latter  consideration  is  one  of 
the  more  productive  changes  in  reducing  execution  time,  because  it  greatly 
reduces  the  time  needed  to  search  the  product  table. 

In  addition,  SPOC  stores  all  the  appropriate  thermodynamic  data 


internal  to  the  code,  rather  than  on  external  tapes  as  done  in  the  original 
NASA  code.  Propellant  ingredients  stored  in  the  code  are; 

Ingredient 

Ingredient 

Name 

Choice  of  Ingredients 

Binder 

BIND 

HTPB  --  Mixture  consisting 

Of  HTPB  nolymer  and  typical 
cure  agent,  plasticizer  and 
bond  agent 

Fuel 

FUEL 

C  - -  Carbon 

ZR  --  Zirconium 

AL  --  Aluminum 

Oxidizer 

OXA,  OXB 

AP  --  Ammonium  perchlorate 
HMX  --  Cyclotetramethylene- 
tetranitramine 

RDX  --  Cyclotrimethylene- 
trinitramine 

Rate  Catalyst 
(Sblid) 

RCATS 

Fe£0  --  Iron  Oxide 

FCH  --  Ferrocene 

Rate  Catalyst 
(Liquid) 

RCATL 

None  available  at  the  present 

Combustion 

Stabilizer 

STAB 

ZRC  --  Zirconium  carbide 

ALOX  -  -  Aluminum  oxide 

ZR  --  Zirconium 

C  --  Carbon 

The  pertinent  properties  for  HTPB  binder  are* 

Formula:  p(7.  133)  H(ll.  150)  0(0.  135)  N(0.  067) 

Enthalpy  (cal/gram  formula  weight):  -6047.  0 
Density  (gm/cu  cm):  0.903 

Input  Information 

Manipulations  within  the  code  to  adjust,  verify  and  normalize  the 
propellant  formulation  use  generic  nomenclature  so  that  the  user  can  have 
a  choice  of  ingredients  for  any  one  constituent  class.  For  example,  two 
oxidizers  may  be  employed  (identified  as  OXA  and  OXB)  and  either  one  of 
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them  may  be  "named"  aa  AP,  or  RDX,  or  HMX.  Furthermore,  each  of  the 
oxidizers  may  be  indexed  to  deaignate  that  they  have  up  to  three  different 
particle  sizes,  to  use,  for  example,  in  propellant  burn  rate  models  that 
combine  amount  and  size,  Particle  size  can  be  input  for  the  AP  oxidizer 
(DIAAP(I),  1=1,  2,  3)  for  use  in  the  combustion  stability  module. 

Up  to  four  other  ingredients  may  be  defined  by  the  user  to  the  code 
to  provide  additional  propellant  formulation  combinations.  The  user  must 
furnish  the  chemical  formula,  enthalpy,  and  density,  along  with  their  iden¬ 
tifying  name.  Then,  this  name  can  ba  employed,  along  with  those  already 
stored  in  the  code,  to  describe  the  propellant.  Table  3  lists  the  species  which 
can  be  employed  in  '.he  chemical  formula  of  a  new  ingredient. 


TABLE  3 

SPECIES  CONTAINED  IN  THERMOCHEMICAL  ANALYSIS  MODULE 


FORMULA 

TEMP  RANGE  (K) 

PHASE 

AL1 

300. 

6000. 

G 

ALl 

932. 

4000. 

L 

AL1 

300. 

932. 

S 

ALl  CL1 

300. 

6000. 

G 

ALl  CL1  0  1 

300. 

6000. 

G 

ALl  CL2 

300. 

6000. 

G 

ALl  CL3 

300. 

6000. 

G 

ALl  CL3 

466. 

1500. 

L 

ALl  H  1 

300. 

6000. 

G 

ALl  H  1  0  1 

300. 

6000. 

G 

ALl  N  1 

300. 

6000. 

G 

ALl  0  1 

300. 

6000. 

G 

ALl  0  1  H  1 

300. 

6000. 

G 

ALl  0  2 

300. 

6000. 

G 

AL2  0  1 

300. 

6000. 

G 

AL2  0  2 

300. 

600G. 

G 

AL2  0  3 

2327. 

4000. 

L 

AL2  0  3 

300.. 

2327. 

S 

C  1 

300. 

6000. 

G 

C  1 

300. 

6000. 

S 

C  1  CL1  O  1 

300. 

6000. 

G 

C  J.  CL4 

300. 

6000. 

G 

C  1  H  1  N  1 

300. 

6000. 

G 

ClHlNlOl 

300. 

6000. 

G 

C  1  H  1  0  1 

300. 

6000. 

G 

C  1  H  2 

300. 

6000. 

G 

C  1  H  2  O  1 

300. 

6000. 

G 

C  1  H  3 

300. 

6000. 

G 

C  1  H  3  CL1 

300. 

6000. 

G 

C  1  H  4 

300. 

6000. 

G 

C  1  N  1 

300. 

6000. 

G 

C  1  N  1  0  1 

300. 

6000. 

G 

C  1  0  1 

300. 

6000. 

G 

C  1  0  2 

300. 

6000. 

G 

C  1  ZP1 

300. 

3805. 

S 

C  2 

300. 

6000. 

G 

C  2  H  1 

300. 

6C00. 

G 

C  2  H  1  CL1 

300. 

6000. 

G 

C  2  H  2 

300. 

6000. 

G 

C  2  H  4 

300- 

6000. 

G 

C  2  N  2 

300. 

6000. 

G 

CL1 

300. 

6000. 

G 

CL1  FE1 

300. 

6000. 

G 

CL1  H  1 

300. 

6000. 

G 

CL1  H  1  0  1 

300. 

6000. 

G 

CL1  N  1  0  1 

300. 

6000. 

G 

CL1  0  i 

300. 

6000. 

G 

CL1  ZR1 

300. 

6000. 

G 

CL2 

300. 

6000. 

G 

CL2  FEl 

300. 

6000. 

G 
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Table  3 


FORMULA 

TEMP 

RANGE  (K) 

PHASE 

CL2  FE1 

950. 

3000. 

L 

CL2  FE1 

300. 

950. 

S 

CL2  ZR1 

300. 

6000. 

G 

CL2  ZR1 

1000. 

2000. 

L 

CL2  ZR1 

300. 

1000. 

S 

CL3  FE1 

300. 

6000. 

G 

CL3  FE1 

577. 

1500. 

L 

CL3  ZR1 

300. 

6000. 

G 

CL3  ZR1 

300. 

2000. 

S 

CL4  FE2 

300. 

6000. 

G 

CL4  ZRl 

300. 

6000. 

G 

FE1 

300. 

6000. 

G 

FE1 

1809. 

4500. 

L 

FE1 

300. 

1809. 

S 

FE1  H  2  0  2 

300. 

6000. 

G 

FE1  H  3  0  3 

300. 

1500. 

S 

FE1  0  1 

300. 

6000. 

G 

FE1  0  1 

1650. 

5000. 

L 

FE1  0  1 

300. 

1650. 

S 

H  1 

300. 

6000. 

G 

H  1  N  1 

300. 

6000. 

G 

H  1  N  1  0  1 

300. 

6000. 

G 

H  1  N  1  0  2 

300. 

6000. 

G 

H  1  N  1  0  3 

300. 

6000. 

G 

H  1  0  1 

300. 

6000. 

G 

H  1  0  2 

300. 

6000. 

G 

H  1  ZRl 

300. 

6000. 

G 

H  2 

300. 

6000. 

G 

H  2  N  1 

300. 

60C0. 

G 

H  2  0  1 

300. 

6000, 

G 

H  2  0  1 

273. 

373. 

L 

H  3  N  1 

300. 

6000. 

G 

N  1 

300. 

6000. 

G 

N  1  0  1 

300. 

6000. 

G 

N  1  O  2 

300. 

6000. 

G 

N  1  ZRl 

300. 

6000. 

G 

N  1  ZRl 

3225. 

6000. 

L 

N  1  ZRl 

300. 

3225. 

S 

N  2 

300. 

6000. 

G 

N  2  0  1 

300. 

6000. 

G 

N  2  O  3 

300. 

6000. 

G 

0  1 

300. 

6000. 

G 

0  1  ZRl 

300. 

6000. 

G 

0  2 

300. 

6000. 

G 

O  2  ZRl 

300. 

6000. 

G 

0  2  ZRl 

2950. 

6000. 

L 

0  2  ZRl 

300. 

2950. 

S 

ZRl 

300. 

6000. 

G 

ZRl 

2125. 

5500. 

L 

ZRl 

300. 

1500. 

S 
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IMPULSE  EFFICIENCY 


The  rocket  motor  impulse  efficiency  is  a  fixed  user  input,  deter¬ 
mined  from  the  empiricisms  dsfined  for  the  SPP(Reference  12) computer 
code,  or  calculated  from  a  user- supplied  model. 

As  employed  in  the  SPOC  ballistic  simulation  module,  impulse 
efficiency  is  the  ratio  of  delivered  vacuum  specific  impulse  at  motor  pressure 
and  expansion  ratio  to  theoretical  specific  impulse  at  identical  conditions. 
Divergence  losses  are  not  included  in  the  efficiency  factor,  but  are  accounted 
for  through  the  thrust  coefficient  calculated  by  the  code.  In  actual  practice, 
the  impulse  efficiency  is  treated  as  a  "thrust  efficiency"  during  ballistic 
simulation.  However,  the  impulse  efficiency  predicted  in  SPP  already 
includes  a  divergence  loss  term.  Therefore,  when  the  SPP  model  is  selected 
by  the  user  (SPPETA=T),  the  divergence  loss  calculation  in  the  ballistic 
module  is  by-passed  by  internally  setting  the  effective  nozzle  exit  half-angle 
to  zero. 

Impulse  efficiency  calculated  by  the  SPP  empiricisms  is  as  follows 
(where  code  nomenclature  is  also  shown) 

ETAISP  =  (ETACSR)(ETACF) 

where 

ETACF  =  1  -  (ETABL)(ETADIf)(ETAKIN)(ETASUB)(ETATP)/100 


and 

ETATP 

ETADl'V 

ETAK1N 

ETABL 

ETACSR 

ETASUB 


is  impulse  loss  effect,  in  %, 

is  impulse  loss  effect,  in  %, 

is  impulse  loss  effect,  in  %, 

is  impulse  loss  effect,  in  %, 
in  the  nozzle 


due  to  two- phase  flow 

due  to  nozzle  divergence 

due  to  finite  rate  reaction  kinetics 

due  to  boundary  layer  buildup 


is  impulse  loss  effect,  in  %,  due  to  C*  efficiency  (i.  e.  , 
same  as  combustion  efficiency) 

is  impulse  loss  effect,  in  %,  due  to  nozzle  submergence. 
Always  equal  to  one  because  submerged  nozzles  are  not  pro¬ 
vided  in  SPOC. 


The  various  losses  are  defined  by  empiricisms  as  described  below. 
ETATP 


The  two-phaBe  flow  loss  is  given  by  the  empiricism 
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ETATP 


where 

£  is  mol  fraction  of  condensed  species  (at  average  motor 

chamber  pressure)  expressed  in  mole  of  condensed 
species  per  100  grams  of  mixture 

D  is  condensed  specie  particle  size  in  microns  and  is  given  by 

.P 

Dp=  0.454  P1/3  {  1/3  [l.e-  004L‘]  [1  +-045  Dt] 


P  is  average  motor  chamber  pressure  in  psia  (PBAR) 

(  iii  nozzle  expansion  ratio  at  ignition  conditions  (ERI) 

R 

is  nozzle  throat  diameter  in  inches  at  ignition  (DTI) 
L*  is  motor  L*  in  inches  at  Ignition  (LSTRI) 

e  is  the  Naperian  base  2.  71828, . . 


The  coefficients  C„,  C,,  Cr,  and  C,  are  determined  as  follows; 
3  4  5  6 


if  Dfc<l;  C3  =  30,  C5  =  1,  C6  =  1 

if  1  <  D.  <  2;  C,  =  30,  Cc  =  1,  C,  =  0.  8 

—  t  —  3  D  o 

if  D,  >  2  &  <  4;  C,  =  44.  6,  C.  =  0.  8,  C ,  =  0.  8 

t  P  3  5  o 

if  D.  >  2  D  <  8;  C  =  34,  C  =  0.  8,  C,  =  0.  4 
t  P  ~  3  5  0 

if  D.  >  2  &  D_  >  8;  C  =  25.  2,  Cc  =  0.  8,  C,  =  0.  33 
t  P  3  5  o 


if£  >  0.09 

C.  =  °.  5 

4 

if  Dfc<l;  C3  =  9,  C5  =  1,  C6  =  1 

if  1  <  D.  <  2;  C  =  9#  C_  =  1,  C,  =  0.  8 
—  t  —  3  5  b 

if  D  >2  &  D  <4;  C.  =  13.4,  C.  =  0.  8,  C,  =  0.  8 
t  P  3  5  b 

if  Dfc  >  2  &  4  <  Dp;  C3  =  1 0.  2,  Cg  =  0.  8,  C&  =  0.  4 

if  D,  >  2  &  D_  >  8;  C  =  7.  58,  C  =  0.  8,  C,  =  0.  33 
t  P  3  5  o 


ETADIV 


The  nozzle  divergence  loss  is  given  by 

etadiv  =  so  [i  -  cosj£J^2L)j 


aEy  =  ALFAEX 


“EX  " 
ALFAEX 


=  ALFA 


Conical  Nozzle 


Contoured  Nozzle 


ETAKIN 


The  reaction  kinetics  loss  is  given  by 

I 


ETAKIN 


f  a 


where 


f 


I  is  theoretical  vacuum  Iap  computed  at  ignition  nozzle 

*V  geometry  and  assuming  Frozen  equilibrium  thermo¬ 

chemistry  (IVACF,  FISP) 

I  is  theoretical  vacuum  Iap  computed  at  ignition  nozzle 

geometry  and  assuming  shifting  equilibrium  thermo- 
chemistry  (IVAC,  SISP) 

if  P  <  200,  A  =  1 

200 

if  P>  200,  A  =  — 

P  is  average  motor  chamber  pressure  in  psia  (PBAR) 

ETABL 

The  boundary  layer  loss  is  given  by 


where 

Cj  is  0.  00365 

C2  is  0.00093  7 

P  ’ "  average  motor  chamber  pressure  in  psia  (PBAR) 

Dfc  is  nozzle  throat  diameter  in  inches  at  ignition  (DTI) 

e  is  the  Naperian  base  2.  71828  . . . 

t  it  tor  burn  time  in  seconds  (TB) 

(  is  nozzle  expansion  ratio  at  ignition  conditions  (ERI) 

ETACSR 

The  c*  efficiency  loss  is  given  by 


ETACSR  = 


K  +  (~y)  (100-K) 


b  (c) 


if  a  <  10,  d  =  a 
if  a  >  10,  d=10 

where 

a  is  the  weight  percent  of  aluminum  in  the  propellant  formulation 
b  =  1.  0 
c  =  1.  003 


The  parameter  K  is 

a  function  of  burning  rate  and  is 

determined  by 

the  following  table. 

Burning  (RBBAR)  is  evaluated  at  motor 

average  pressure 

Burning  Rate 

Burning  Rate 

(in/ sec) 

K 

(in/sec) 

K 

<0.11 

91.4 

0.50 

98.  6 

0.  11 

91.4 

0.  60 

98.  9 

0. 12 

9^.1 

0.  70 

99. 1 

0.13 

94.  0 

0.  80 

99.2 

0.  14 

94.6 

0.90 

99.3 

0.  15 

95.  1 

1.  00 

99.4 

0.  16 

95.6 

1.20 

99.  6 

0.  17 

96.  0 

1.40 

99.7 

0.  18 

96.4 

1.  60 

99.  8 

0.  19 

96.  7 

1. 80 

99.9 

0.20 

97.  0 

2.  00 

100.  0 

0.30 

97.  7 

>2.  00 

100.  0 

0.40 

98.2 
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Values  for  the  parameters  used  to  calculate  the  various  loss  factors 
are  generated  at  numerous  places  in  the  code.  Table  4  lists  the  param¬ 
eters  (according  to  both  the  calling  and  called  arguments)  and  their  sources. 
The  primary  influence  was  the  need  to  use  data  from  the  just  completed 
evaluation  (the  precedihg  pass  through  COMP)  in  order  to  not  have  an  itera¬ 
tion  within  one  COMP  evaluation.  For  example,  the  proper  cha  nber  pres¬ 
sure  to  use  in  estimating  ETAISP  is  the  one  at  which  the  motor  will  operate; 
but  since  pressure  and  thrust  simulations  are  performed  concurrently,  the 
average  pressure  is  not  available  when  the  estimate  for  ETAISP  is  needed. 
Thus  the  technique  of  using  pressure  from  the  just  completed  simulation  was 
adopted.  The  error  usually  should  be  small  (at  least  not  excessive),  and 
as  the  optimum  design  is  approached  the  small  step-sise  by  PATSH  will 
produce  only  small  changes  in  basic  parameters.  Therefore  the  conditions 
for  the  just  completed  evaluation  should  be  essentially  identical  to  the  cur¬ 
rent  evaluation.  As  seen  in  Table  4,  there  are  other  parameters  treated 
in  the  same  manner;  unless  otherwise  stated,  the  data  come  from  the  cur¬ 
rent  evaluation. 


TABLE  4 


SOURCES  OF  DATA  USED  IN  SPP  PREDICTION 
OF  IMPULSE  EFFICIENCY 


Calling 

Argument 

Called 

Argument 

Definition  and  Units 

TB 

TB 

Burn  time  (sec)  from  zero  to  time^at  which 

99.  5%  of  the  propellant  has  been  consumed. 

Taken  from  preceding  ballistic  simulation  at 
low  temperature  (>f  a  two-temperature  problem), 
otherwise  from  the  single  temperature  simulation. 

One  second  on  initial  pass. 

P 

PBAR 

Average  pressure  (psia)  over  TB,  taken  at  same 
conditions  as  TB,  from  preceding  pass  through 
COMP.  Input  PC  (default  to  1000  psia)  on  initial  pass 

ALFA 

ALFA 

Half-angle  (deg)  at  entrance  to  exit  cone  in 
contoured  nozzle.  Equal  to  ALFAEX  in  conical 
exit  cone. 

ALFAEX 

ALFAEX 

Half-angle  (deg)  at  exit  of  exit  cone. 

NOZER 

ERI 

Nozzle  initial  expansion  ratio. 

RATE 

RBBAR 

Burning  rate  (in/sec)  at  pressure  PBAR. 

Burn  rate  subroutine  called  with  pressure  P 
to  obtain  this  value. 

LSTRI 

LSTRI 

Initial  L*  (in).  Ratio  of  initial  port  volume 
from  preceding  evaluation  to  initial  throat 
area.  Port  volume  includes  only  that  part  of 
chamber  occupied  by  propellant.  Estimate 
for  first  pass  through  COMP  is 

00  (RMOT  OS)2  ( LMOT  MX )  /  3 . 

DTI 

DTI 

Initial  throat  diameter  (in). 

IVACF 

FISP 

Theoretical  specific  impulse  (lbf-sec/lbm)  at 
frozen  equilibrium,  current  nozzle  expansion 
ratio,  preceding  chamber  pressure;  from  sub¬ 
routine  TCHEM. 

IVAC 

SISP 

Theoretical  specific  impulse  (lbf-sec/lbm) 

at  shifting  equilibrium,  current  nozzle  expan¬ 
sion  ratio;  at  preceding  chamber  pressure; 
from  subroutine  TCHEM. 
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Table  4  (Continued) 

SOURCES  OF  DATA  USED  IN  SPP  PREDICTION 
OF  IMPULSE  EFFICIENCY 

% 

Calling  Called  *  » 

Argument  Argument _ Definition  and  Unita  >, _ 

MOLCND  MOLFR  Mole  fraction  of  condensed  species  (*noles 

per  100  gm  of  mixture).  From  subroutine 
TGHEM  at  preceding  chamber  pressure. 

PCAL  PCAL  Percent  aluminum  (%)  in  propellant. 


NOZZLE  STRUCTURAL  AN D  ABLATIVE  THICKNESS 


Each  machine  access  will  require  the  user  to  specify  a  nozzle  type 
to  be  considered.  The  user  may  choose  one  from  among  six  nozzle  types 
shown  by  Figures  20  through  25. 

All  nozzle  types  consist  of  an  inner  layer  of  insulating  material 
supported  by  an  outer  structural  member.  For  Nozzle  Types  1  and  2  and 
the  exit  cone  of  Type  6,  the  thickness  of  the  insulating  material  is  established 
by  the  erosion,  char  and  thermal  penetration  depths  calculated  in  the  code. 

For  Nozzle  Types  3,  4,  5  and  the  entrance /blast  tube  of  Type  6,  the  thictuiess 
of  the  insulating  material  is  established  by  user  inputs  that  define  the  outer 
contour  of  the  structural  material,  which  in  turn  defines  the  outer  contour  of 
the  insulating  material.  Calculations  of  erosion,  char  and  thermal  penetration 
depths  are  still  made  for  the  latter  nozfelee,  and  the  results  are  used  to 
determine  if  sufficient  insulation  material  has  been  provided. 

Thickness  requirements  of  the  structural  members  are  calculated  in 
the  code  for  all  the  nozzles  except  for  Type  3  and  the  supersonic  portion  of 
Type  4. 

Up  to  three  different  insulating  materials  can  be  specified  for  Nozzle 
Types  1,  2  and  6.  The  boundary  between  materials  is  defined  by  a  user- input 
area  ratio  (ARl  in  the  entrance  section,  AR2  in  the  exit  section).  By  definition. 
Nozzle  Type  3  and  4  have  only  one  insulating  material.  Again,  by  definition, 
the  boundary  between  insulating  materials  No,  1  and  No.  2  occurs  at  the  aft 
end  of  the  blast  tube  on  Nozzle  Type  5  and  Type  6.  The  boundary  between 
insulating  material  No.  2  and  No.  3  in  Ntozzle  Types  1,  2  and  6  will  occur 
at  the  AR2  input  by  the  user;  however  the  conical  ramp  supporting  the  insert 
of  Nozkle  Type  2  will  be  positioned  to  connect  the  two  boundaries  regardless 
of  the  relative  magnitudes  of  ARl  and  AR2. 

Two  structural  support  materials  can  be  specified.  Material  No.  1 
is  for  the  entrance  and  throat  regions;  No.  2  is  for  the  exit  region.  The 
boundary  between  Structural  Materials  No.  1  and  2  occurs  a  distance  XSTRAN 
downstream  of  the  boundary  between  Insulating  Materials  No.  2  and  3. 

The  stagnation  pressure  from  which  the  local  static  pressure  is 
calculated  is  MEOP  (maximum  expected  operating  pressure)  determined 
from  the  high  temperature  ballistic  simulation;  thus  all  pressure-dependent 
analyses  in  the  nozzle  subroutine  Include  an  inherent  degree  of  conserva¬ 
tism  by  the  use  of  the  MEOP,  which  is  the  upper  three- sigma  maximum 
pressure.  Additional  conservatism  is  included  in  motors  with  either  a  pro¬ 
gressive  or  regressive  pressure  history  which  have  an  average  pressure 
less  than  the  maximum  on  which  MEOP  is  based. 

Burn  time  used  in  the  nozzle  analyses  is  the  nominal  burn  time  deter¬ 
mined  from  the  low  temperature  ballistic  simulation. 
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Establishing  Internal  Contour 

The  program  uses  the  parameters  shown  in  Figure  6 6  to  describe  the 
internal  contour  of  the  nossle  In  terms  of  an  X«R  coordinate  system  with  its 
origin  at  the  nosale  entrance. 


■Q^  =  ALFA 


Tangent  Point 


a2  =  ALFAEX 


Figure  67.  Condition*  for  Contoured  Expansion  Section 

The  parameters  shown  in  Figure  67  are  used  to  determine  whether 
the  nozzle  is  to  have  either  a  conical,  elliptic,  parabolic  or  hyperbolic  exit 
cone.  This  decision  is  made  in  the  following  way. 


Condition 


Ofj  =  a2  (input) 


r*  ~t2  r  *— 1  2 

T»n°jl  K 
“2J  >  L\ 

p-  aif  r  \ 

L™  “0  '  L\ 

P"  M  2  .  Pe 

[r»n«2J 


Decision 


Conical  Exit  Section 


Elliptic  Exit  Section 


Parabolic  Exit  Section 


Hyperbolic  Exit  Section 


The  inner  contour  is  described  by  an  X- R  array,  where  the  X- 
coordinates  are  established  by  a  user-input  incremental  X.  The  procedure 
to  establish  the  inner  contour  is  the  same  for  all  six  noaxle  types.  The  only 
difference  comes  in  defining  the  blast  tube  length  for  Types  5  and  6;  there 
a  radius  equal  to  DTI  ^/aRi  is  held  constant  until  a  length  LBT  is  reached. 
Figure  66  shows  general  contour  information. 


Establishing  Insulation  Thickness 


After  describing  the  internal  contour,  the  insulation  thickness  is 
determined  for  each  point  in  the  nozzle  by  considering  the  amount  needed  for 
erosion,  char,  thermal  protection  and  safety  factor. 


The  thickness  of  insulation  needed  due  to  erosion  is  determined  at 
each  X-R  coordinate  from  a  mathematical  model 

Ca  S 

r  =  C.  P  4  M  (1  +  C,  sin  Otl  +  C,P  +  C,  ( 1 ) 

e  1  2  3  6 

%  =  <V  (,b>  (2) 

where 


rg  =  local  erosion  rate  (in/sec) 

T  =  local  erosion  thickness  (in) 

e 

t  =  motor  burn  time  (sec) 

b 

P  -  local  static  gas  pressure  (psia) 

M  =  local  gas  Mach  number 

a  =  local  angle  between  nozzle  internal  contour  and 

nozzle  centerline  (deg) 

Cj,  C2 - =  parameters  that  are  a  function  of  the  insulating 

material  and  its  location  in  the  nozzle 


Values  for  P  and  M  are  calculated  from  the  isentropic  flow  relations. 

Local  values  for  ft  come  from  the  contour  routines.  Cj,  C£ - are 

determined  from  a  statistical  curve  fit  of  erosion  data  for  the  material  under 
consideration,  input  by  user. 


The  thickness  of  insulation  needed  due  to  char  is  determined  at  each 
X-R  coordinate  from  the  following  relation  found  in  Reference  13  for  an  ablating 
surface 


Nozzle  Geometry  Showing  X-R  Coordinates  Described  in  Output  of  Nozzle 
Analysis  Routine 


x  is  limited  to  a  maximum  of: 


r  =  0.  79  ^ 
c 


a  t 
c  b 


-In 


1- 


T  -T 
char  amb 


-  +0.574 


T  -T  , 
vapor  amb 


(4) 


where 


char  thickness  (in) 


a 

c 

k 

c 

T 

char 

T 

vapor 

T 

amb 

h 


char  thermal  diffusivity  (sq  in/sec) 

char  thermal  conductivity  (BTU/in- sec- °F) 

material  char  tamparatura  (*F) 

material  vaporization  temperature  (#F) 

ambier^oyL||itiayo|emperatur  e  of  material  (*F). 

surface  heat  transfer  coefficient  (h— for 
extreme  conservatism)  (BTU/sq  in-sec-°F) 


The  limiting  value  is  obtained  from  a  curve  fit  of  the  Hottel  chart  for  the 
midplane  temperature  of  a  large  slab  (found  in  Reference  14). 

The  thickness  of  insulation  required  as  a  thermal  barrier  for  the 
nozzle  structural  material  is  determined  from  the  following  relationship 


rb  =  - 


r°b ‘b  “ 

r  t  -  t 

“  alow  amb 

In  1  ~  t'£  -  X 

|  v*  po  r  amb 

1 

*-u 

+ 

_ 1 

C,  +  (C,  -  1. 0)  T  (5) 
b  b  e 


where  Tr  if  limited  to  a  maximum  of 
b 


fb  =  cb 


* 

,  2 

'  a 

t. 

kb 

0.79  « 

c 

b 

+  0  574  — 

- 

T  -T 

alow  amb 

L  h  J 

l 

“  “  T  -T 

.char  amb 
h  1  J 

t 

(b) 


“fti 


+  (Cb'  ,-0,T. 
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and  where 


b 

CL 


‘alow 


thermal  barrier  thickness  (in) 
barrier  thermal  diffusivity  (sq  in/sec) 
barrier  thermal  conductivity  (BTU/in-sec-°F) 
thermal  barrier  safety  factor  (provided  by  user) 
allowable  temperature  (®F)  of  structural  material 


The  limiting  value  is  also  obtained  from  a  curve  fit  of  the  Hottel  chart  in 
Reference  14. 

The  total  insulation  thickness  is 


r.  =  r  +  r  +  r. 

i  e  c  b 


(7) 


noting  that  a  safety  factor  is  included  in  the  thermal  barrier  increment. 

Separate  data  sets  of  input  parameters  must  be  furnished  for  each 
of  the  insulating  materials  inherent  with  a  particular  nozzle  type  (three  for 
Type  1,  one  for  Type  3,  one  for  Type  4,  etc.  ). 


Cl*  C2*  * '  C6 

Cb 

T 

vapor 

Tchar 


The  thickness  of  the  insulation  in  the  entrance  section  of  Nozzle 
Types  4,  5  and  6  is  input  by  the  user  (TENT);  however,  erosion,  char  and 
thermal  barrier  thickness  requirements  are  also  calculated.  The  final 
thickness  is  the  greater  of  the  two  (input  or  calculated). 

Establishing  Structure  Thickness 

Nozzle  Type  1  and  Type  2 

Structural  requirements  for  nozzle  Type  1,  Type  2  and  the  exit  cone 
of  Type  6  are  found  with  the  following  analysis.  The  longitudinal  and  radial 
forces  the  nozzle  must  support  at  any  point  are  calculated.  (See  Figure  69.  ) 
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Fx  ia  the  longitudinal  force  found  by  summing  the  pressure  forces  from  the 
exit  plane,  Fy  is  the  hoop  load  caused  by  the  local  static  pressure. 


Figure  69.  Forces  Acting  on  Nozzle 


The  structural  thickness  needed  to  withstand  these  forces  is  cal¬ 
culated  from  the  following  relations.  Load  carrying  capability  of  the 
insulating  material  is  ignored. 

(1)  Longitudinal  Stress 

F 

f  -  - 5 - 

1  2  ir  r  cos  a  fl 

c 


where 


Tj  =  thickness  (in)  required  by  longitudinal  stress 

=  longitudinal  force  (lbf)  at  individual  station, 
summed  from  exit  plane 

R  =  local  radius  (in) 


\ 

I 

i  \ 

\  \ 


■where 


where 


where 


a 


'2 

E 


'4 

P 

fi 


t 


angle  (deg)  between  nozzle  surface  and  nozzle 
centerline 

compressive  yield  strength  (pai)  of  structural 
material 


(2)  Buckling  Stress 


V5  Fa 
2  it  cot 


cos  a  E 


(9) 


thickness  (in)  required  by  buckling  stress 

modulus  of  elasticity  (psi)  of  the  structural 
material 


(3)  Shear  Stress 


F  sin  a 

f  -  -2L - 

3  2  It  R 


(10) 


thickness  (in)  required  by  shear  stress 

shear  yield  strength  (psi)  of  structural 
material 


(4)  Hoop  Stress 


r  =  ££• 

4 


(ii) 


thickness  (in)  required  by  hoop  stress 
local  gas  pressure  (psia) 

tensile  yield  strength  (psi)  of  structural  material 
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The  thickness  of  the  structural  material  is  the  maximum  of  fj,  r^, 
fy  or  r  .  The  input  structural  safety  factor  is  then  applied  to  the  maximum 
thickness. 

Nozzle  Type  3 

A  structural  analysis  is  not  performed  on  Type  3  nozzles  because  the 
ablative  adds  significantly  to  the  load-carrying  capability  of  the  assembly. 
Generally  speaking,  the  capacity  of  the  structure  is  not  a  critical  item. 

When  a  Type  3  nozzle  is  employed  with  an  Aft  Closure  Type  2,  the  structure 
thickness  is  made  equal  to  the  case  thickness  (TCASE).  When  it  is  com¬ 
bined  with  an  Aft  Closure  Type  1,  the  support  structure  thickness  is  a  user 
input  TSTR3. 

Nozzle  Type  4 

A  structural  analysis  is  not  performed  on  the  reduced-diameter  aft 
section  of  Type  4  nozzles;  its  thickness  is  a  user  input  TSTR4.  Thickness 
of  the  structure  that  forms  the  entrance  section,  TENTS,  (Figure  70)  is 
estimated  through  the  thin-wall  pressure  vessel  relationship 


where 


TENTS 

RENT 


(P)(RENT) 

(FTY)(FSTRUS) 

RNl/sin  (90-ALFAEN) 


(12) 

(13) 


P  = 
RENT  = 
FTY 

FSTRUS  = 


MEOP  (psia) 

Effective  radius  (in);  See  Figure  70 
Structural  material  No.  1  tensile  yield  strength 

(pal) 

Structural  material  No.  1  safety  factor 


Nozfcle  Type  5  and  6 

Thickness  of  the  structure  that  forms  the  entrance  to  the  blast  tube 
is  calculated  with  Eqs  (12)  and  (13).  Structure  in  the  blast  tube  is  calculated 
with 


where 


T  RTS  -  P(KUP)  *  R3  (KUP) 

FTY  *  FSTRUS 

TBTS  =  Thickness  (in)  of  blast  tube  structure 

P(KUP)  =  Static  pressure  (psia)  at  entrance  to  throat 
insert  (Station  KUP) 

R3(KUP)  =  Outside  radius  (in)  of  insulation  material  along 
blast  tube  (Station  KUP) 
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Figure  70.  Basis  for  Structural  Analysis  of  Entrance  Section 

of  Nozzle  Types  4,  5  and  6 
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Structure  thickness  for  the  exit  cone  of  Nozzle  Type  6  is  determined 
with  Eqs  (8)  -  (11),  as  were  Types  1  and  2. 

Interfaces  with  Motor  and  Geometric  Verification 


The  outermost  contour  of  Nozzle  Type  1,  Type  2  and  the  exit  cone 
of  Type  6  are  a  result  of  separate  calculations  for  erosion  depth,  char  depth, 
thermal  barrier  and  structural  thickness,  all  performed  at  a  number  of  node9 
(located  by  X-R  coordinates)  along  the  internal  surface  of  the  motor.  Cal¬ 
culations  at  each  node  are  independent  of  one  another,  and  so  there  is  some 
"waviness"  in  the  outer  contour  (except  for  behind  the  throat  of  Type  2  nozzle). 
It  is  recognized  that  a  nozzle  would  not  actually  be  built  with  this  contour,  but 
the  purpose  of  these  analyses  is  to  obtain  a  reasonably  accurate  estimate  of 
the  nozzle  size,  and  to  that  extent,  the  analyses  are  appropriate. 

As  described  above,  the  outer  contours  of  Nozzles  Type  3,  Type  4, 
Type  5  and  the  blast  tube  of  Type  6  are  established  by  various  user  inputs. 

How  the  interface  between  the  user-established  outer  contour  and  the 
analysis-established  inner  contour  is  controlled  is  described  in  later  para¬ 
graphs. 


The  primary  interface  between  the  motor  and  nozzle  is  at  the  aft  case 
opening  where  the  case  opening  radius  RNOZEN  must  mate  with  the  nozzle 
entrance  radius  RNl.  For  any  given  problem,  a  value  of  RNOZEN  is  cal¬ 
culated  at  every  design  evaluation  (every  pass  through  COMP),  whatever 
the  combination  of  propellant  grain  configuration  and  aft  closure  type  might 
be.  The  definition  of  the  various  possibilities  of  RNOZEN  are  illustrated  in 
the  figures  included  as  part  of  the  grain  cojafigu ration  discussion. 

A  series  of  geometric  validations  are  made  prior  to  ballistic  simula¬ 
tion  to  achieve  the  proper  motor/nozzle  interface  and  to  assure  compatibility 
of  other  nozzle  dimensions  (Figures  71  ,  72  and  73).  These  comparison 

guarantee; 

(1)  Radius  of  boundary  between  Insulation  Materials  No.  1  and  No.  2 
(RT ARl)  is  less  than  RNOZEN  (Nozzle  Types  1.  2,  5  and  6).  Note  that 
for  Types  5  and  6,  this  material  boundary  radius  corresponds  to  the  inside 
radius  of  the  blast  tube. 

(2)  Exit  radius  (RE)  is  greater  than  the  radius  of  the  boundary _ 

between  Insulation  Material  No.  2  and  Insulation  Material  No.  3  (RT AR2) 
for  Nozzle  Types  1,  2  and  6. 
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Geometric  Verification  of  Nozzle  Types  3  and  4 


write  "RT>/nrr*  rznmx1 


(3)  Throat  radius  (RT)  is  less  than  aft  case  opening  (RNOZEN)  for 
Nozzle  Types  3  and  4. 

(4)  Exit  radius  (RE)  is  greater  than  throat  radius  (RT)  for  Nozzle 
Types  3  and  4. 

(5)  Exit  radius  (RE)  is  less  than  the  nozzle  entrance  radius  (RNOZEN) 
for  Nozzle  Type  3. 

(6)  Exit  radius  (RE)  is  less  than  the  inside  radius  of  the  blast  tube 
support  structure  (RBTO-TSTR4)  for  Nozzle  Type  4. 

(7)  Nozzle  entrance  radius  (RN1)  is  equal  to  the  aft  case  opening 
(RNOZEN)  for  Nozzle  Types  1,  2,  5  and  6. 

(8)  Exit  half-angle  (ALFAEX)  is  less  than  the  expansion  section 
entrance  angle  (ALFA)  for  contoured  nozzles. 

When  the  ballistic  simulation  is  completed,  the  nozzle  thermal  and 
structural  analyses  are  performed.  These  results  are  used  to  determine 
dimensional  compatibility  between  insulation  and  support  structure  in 
Nozzle  Types  3,  4,  5  and  6.  Figure  74  illustrates  the  check  made  at  the 
exit  plane  of  Nozzle  Types  3,  4  and  5;  a  penalty,  OBEXIM,  is  calculated  if 
the  margin  EXINSM  is  less  than  zero.  For  the  blast  tube  of  Nozzle  Types 
5  and  6,  the  outside  etructure  radius  is  found  by  summing  the  required 
insulation  thickness  (erosion,  char,  thermal  barrier)  and  the  required 
support  structure  thickness  with  the  inside  -adius;  if  this  total  is  greater 
than  the  user-input  RBT05  (or  RBTOb),  a  penalty  is  calculated  (OBBTOB 
or  0BBT06). 

The  outside  radius  of  the  exit  section  (RAO)  is  also  determined  for 
Nozzle  Types  1,  2  and  6  after  the  nozzle  analyses  are  performed.  If  RAO 
is  greater  than  a  user  limit,  a  penalty  is  calculated  (OBJDEO).  The  outside 
exit  radius  limit  is  input  as  the  ratio  of  nozzle  diameter  to  motor  diameter 
(NTMR).  Thus  tne  nozzle  exit  size  can  be  adjusted  in  concert  with  the 
motor  diameter  if  the  latter  is  one  of  the  adjusted  parameters  in  an  optimi¬ 
zation  problem. 

Another  check  is  on  the  length  of  the  blast  tube  (i,  e.  .  the  reduced 
diameter  aft  section)  of  a  Type  4  nozzle.  If  the  length  calculated  in  the 
analysis  (LBT4)  is  not  equal  to  the  required  length  LBT4RQ,  a  penalty  is 
calculated  (OBLBT4). 


PRESSURE  VESSEL  STRUCTURAL  ANALYSIS 


Thicknesses  of  structural  material  required  to  withstand  maximum 
expected  operating  pressure  (MEOP),  plus  a  safety  factor,  are  calculated 
with  conventional  relationships  (Reference  IS).  To  account  for  the  several 
case/closure  arrangements,  the  bfific  stress  relationships  are  employed 
for: 

o  Flat  plate  unrestrained  on  outer  edge,  such  as  when 
closure  is  held  in  place  by  retaining  ring  (Forward 
Closure  Type  2) 

o  Flat  plate  restrained  on  outer  edge,  such  as  when  closure 
is  integral  with  tubular  portion  of  case  (Forward  Closure 
Type  3) 

o  Elliptical  domes  for  both  forward  and  aft  closures  (Type  1) 
These  calculations  are  performed  in  the  subroutine  CASEAN. 

Two  structural  design  pressures  are  calculated  immediately  after 
the  ballistic  simulation;  one  is  for  use  v/ith  yield  tensile  properties  and  the 
other  is  for  use  with  ultimate  tensile  properties.. 

PYIELD  =  (FSYLD)(MEOP)  (1) 

PULT  =  (F3ULTKMEOP)  (2) 


where  FSYLD  =  Factor  of  safety  fer  yield  conditions 

FSULT  =  Factor  of  safety  for  ultimate  conditions 
MECP  -  Maximum  expected  operating  pressure  (psia); 

upper  three- sigma  maximum  pressure  at  high 
temperature  firing 


A  test  is  made  to  determine  whether  the  yield  condition  or  the  ultimate 
condition  is  the  more  critical;  th«J  decision  depends  on  the  relationship  between 
the  two  safety  factors  and  the  material  yield  and  ultimate  strengths 


F 


(FSYLD)(FTUC) 


(3) 


where  Fl'YC  =  Case  structural  material  tensile  yield  strength  (pai) 

FTUC  =  Case  structural  material  tensile  ultimate  strength  (psi) 


If  F>1,  the  ultimate  condition  is  more  critical  (i.  e.  ,  when  designing  to 
ultimate  condition  such  that  the  ultimate  factor  of  safety  is  FSULT,  the 
resultant  yield  factor  of  safety  will  bn  greater  than  that  required,  FSYLD). 
If  F<  1,  the  reverse  situation  prevails.  Note  that  this  test  uses  case 
material  properties  and  the  decision  is  applied  even  to  the  forward  closure 
that  is  a  separate  part  (Type  2). 
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Ellipsoidal  Dome  Closures,  Forward  and  Aft  (Type  1) 


Thickness  of  an  ellipsoidal  closure  is  calculated  from  Case  5,  Table 
XIII  (Reference  15),  considering  conditions  at  the  centerline  (Figure  75). 


The  radius  R.  can  be  shown  to  be 

^  1 

R2  =  a  [  (x/a)2  (1-K2)  +  K2  J  (1) 

2 

Where  K  =  a/b.  When  x  =  zero,  =  a  /b  =  aK.  Thus  in  terminology  used 
in  the  code,  R.,  becomes 

4. 

R  =  (RCI)(BETA2F)  (2) 

where  RCI  =  case  inside  radius  (in) 

BETA2F  =  ellipse  ratio  of  dome  internal  surface  =  a/b 


'  ■",V”'  V  n TfCTTr-.V^y  VY 


From  the  above  quoted  Reference  15  case,  hoop  and  meridional  stress  are 
equal  at  the  center,  so  that  the  required  closure  wall  thickness  (TCLOF  or 
T CLC A)  is 


TCLOF  (or  TCLOA) 


<2)(STRESS) 


where  P  =  critical  pressure  (psia),  ultimate  or  yield,  selected  as 

described  above 

R  =  radius  (in)  from  Eq  (2) 

STRESS  =  case  material  strength  (psi),  ultimate  or  yield, 
selected  as  described  above. 

If  the  thickness  calculated  with  Eq  (3)  is  leas  than  the  thickness  of  the 
case  cylindrical  wall  (TCASE),  it  is  set  equal  to  TCASE.  Mathematically, 
the  condition  of  TCLOF  <  TCASE  will  occur  when  BETA2P<2;  however, 
manufacturing  experience  has  shown  that  it  usually  involves  extra  expense  to 
provide  the  thinner  closure. 

Flat  Plate  Forward  Closure  Not  Integrally  Attached  (Type  2) 


For  a  flat  plate  loaded  as  shown  in  Figure  76,  the  required  thickness 
is  related  to  radial  stress  at  the  center  by  (Case  1,  Table  X,  Ref.  15) 


8ff  m  t 


L  8  t  j 


3m  +  1 


from  which 


/  3  Pr  v  (3m  +  1] 

V  8  F.  " 


=  TCLOF 


Nomenclature  is  in  Table  5. 


Note:  No  moments  exiEt 
at  plate  edges 


Figure  76.  Unrestrained  Flat  Plate  Forward  Closure 

This  type  of  forward  closure  stress  analysis  applies  when  the  flat  plate  is 
keyed  (or  similarly  attached)  to  the  case  wall.  Case  strength  levels  are 
used  to  calculate  this  closure  thickness. 

Flat  Plate  Forward  Closure.  Integral  with  CaBe  | 


Stress  in  this  closure  is  found  by  superimposing  radial  stresses,  Sr, 
as  defined  in  Reference  15  (radial  at  center  of  closure). 


Pressure  (Case  1,  Table  X) 


S  =  - 
r 


3  PR^  (3m  +  I] 
8  m  t, 


Moment  (Case  12,  Table  X) 

s  =  -6-^ 

r  t2 


Shear 

■,-f 

where  M  and  V  are  found  from  Case  24,  Table  XIII. 


Superposition  and  then  combining  terms  results  in 


TABLE  5 


NOMENCLATURE  FOR  PRESSURE  VESSEL 
STRUCTURAL  ANALYSIS 


P:  motor  pressure  in  psia,  yield  or  ultimate  conditions,  PULT,  PYIELD 

r:  plate  radius  in  inches  (same  as  motor  case  inside  radius),  RCI 

m:  reciprocal  of  Poisson's  ratio  for  case  material,  XM 

E:  modulus  of  elasticity  for  case  material,  MODCAS 

t:  forward  closure  thickness  in  inches,  TCLOF 

V:  Poisson's  ratio  for  case  or  closure  material,  PRCAS,  PRCLO 

R:  radius  (in  inches)  from  motor  centerline  to  center  of  case  thickness, 

CAPR 

t  :  case  wall  thickness  in  inches,  TCASE 

c 

M:  moment  applied  to  case  lip  in  inch- pounds  per  inch  of  circumference 

V:  radial  shear  load  applied  at  case  lip  in  pounds  per  inch  of  circumference 

F^:  tensile  strength  of  case  material- --either  ultimate  or  yield,  FTUC, 

FTYC 

R  :  case  outside  radius  in  inches,  RMOTOR 

m 

a;  ellipse  semi-minor  diameter  in  inches,  A1A  or  A1F 

b;  ellipse  semi-major  diameter  in  inches,  BlA  or  BlF 

fr:  radial  stress  in  psi  (tension) 
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where 


<*2  ♦  «3>/<*4  -  Vj 

[<x2+  x3)/(x4-x5)]-  x2 


3  PR2  (3m  +  I) 

*1  =  8mt2 


*  =  12  PE3.'/'  D  (]■/) 

2  4  (1  +  V)  E  t-* 


2  3 

2  PR  \  t  D _ 

tc  (1  -  V/2)  [  Et  +  2RD\3  (l-y)j 

2X+  24R\ZD  (1-v2) 

Et3  (1  +  V) 


(9) 


(10) 


(11) 


(12) 

(13) 


K  Et 


A.  = 


Et  +  2D\  R  (l-V) 
,  -,0.25 

.  3  P~v  l 


2  2 
R  t 

c 


Et; 


D  = 


12  (1-p  ) 


(14) 


(15) 


(16) 


Eq  (9)  is  solved  iteratively  for  t  (i.  e. ,  TCLOF  and  TCLOA)  until  the  calculated 
radial  stress  Sf  is  within  0.  1%  of  the  critical  strength  (either  FTUC  or  FTYC). 

Case  Cylindrical  Section 


The  case  thickness  required  to  withstand  the  predicted  design  pres 
sure  (TCREQ)  is  calculated  by 


TCREQ 


(P)(Rm) 


tx 


(17) 


This  is  compared  with  the  current  value  of  case  thickness  (TCASE)  that  was 
used  to  determine  the  propellant  external  dimensions.  If  TCASEX  TCREQ, 
a  penalty  is  calculated, 
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PROPELLANT  STRUCTURAL  ANALYSIS 


Strain  imposed  on  the  propellant  grain  due  to  low- temperature 
storage  and  low -temperature  ignition  is  calculated  with  "plane- strain" 
relations  (Reference  16)for  the  cylindrically  perforated  cavity.  Strain 
in  the  valley  of  a  star  or  slotted  configuration  is  found  by  application  of 
strain  concentration  factors  (Reference  17)to  strain  calculated  for  equiv¬ 
alent  cylindrically  perforated  ports. 

It  was  decided  to  use  plane  strain  analyses  of  circular  port 
geometries  to  determine  these  strains  and  to  modify  these  with  appropriate 
concentration  factors  when  the  port  is  not  circular.  Such  analyses  should 
be  conservative  for  storage  condition  loading,  in  that  plane  strain  analy¬ 
ses  predict  higher  hoop  strains  (and  stresses)  than  do  three-dimensional 
analyses.  Also,  "lobes"  of  propellant  protruding  into  the  bore  restrict 
hoop  strain  somewhat,  but  this  effect  is  ignored.  For  ignition  conditions, 
plane  strain  analyses  are  expected  to  be  accurate  except  at  transitions 
between  irregular  ports  and  circular  ports.  At  the  transitions  the  actual 
hoop  stress  and  hoop  strain  are  bounded  from  above  by  plane  stress  values 
and  from  below  by  plane  strain  results. 

For  total  hoop  strain  imposed  by  low  temperature  storage  (Table  I, 
Reference  15): 

r 

'9t  =  -  <1  +l/)te)(AT)  [‘—OLjittii - ][  l-2„+(f)2  £,2 


where 

V 


a 


a 

c 


at  = 


a 


b 


r 


Poisson's  ratio  of  propellant  (in/in),  PRP 

Poisson's  ratio  of  case  (in/in),  PRCAS 

Coefficient  of  thermal  expansion  of  propellant  (in/in/ 4F), 
ALPHAP 

Coefficient  of  thermal  expansion  of  case  (in/in/°F), 
ALPHAC 

Conditioned  temperature  minus  strain-free  temperature 
(° F),  DELT  =  TLO-SFTEMP 

Bore  radius  (in),  R2 

Propellant  outBide  radius  (in),  RF 

Radius  at  which  calculations  are  being  made  =  a  when 
bore  strains  are  being  calculated 
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Intermediate  function  defined  in  E'-i  (2)  for  thermal 
loading  of  the  grain,  SIGTHM. 


n. 


(i 


b 


+  i 


[<!•-*][>-/] 


(1  +  V)(h)(Ec) 


(2) 


E 

c 

h 


Equilibrium  modulus  of  prc  pellant  for  use  in  thermal 
loading  (psi),  MO  DPT 

Tensile  modulus  of  caBe  (psi),  MODCAS 
Case  wall  thickness  (in),  TCASE 


Total  strain  found  by  Eq  (1)  is  composed  of  mechanical  (Em)  and  thermal 
(Et)  components;  the  mechanical  component  is  compared  to  strain  endurance 
capability  of  the  propellant  to  determine  the  structural  margin  of  safety 


*  m  '  *et  t 

=  '  0t  -  CV  (AT)  (3) 

where  *  -  Mechanical  hoop  strain  imposed  on  propellant  due 

to  low  temperature  storage  (in/in),  EPT 

€  =  Total  hoop  strain  imposed  on  propellant  -lue  to  low 

temperature  storage  (in/in) 

*  =  Thermal  strain  (in /in) 


For  hoop  strain  imposed  by  ignition  pressurization  at  low  temper¬ 
ature  (Table  II,  Reference  1),  which  are  superimposed  on  thost  strains 
due  to  thermal  shrinkage 


f  _(1  +  i/)(P) 

** 

[l-2p+  (~) 

.  r  ,2  21 

f6P 

r  k2 

L  '  1 

i »  _ 

E 

PP  J 

(4) 


where 


«  =  Hoop  strain  due  to  pressurization  (in/in),  EPP 

Sp 

P  =  Ignition  pressure  (psia),  PIGN 

E  =  Modulus  of  procellant  appropriate  for  ignition  pres- 
^  Bur.lzation  rate  and  temperature  (psi),  MODPP 


Intermediate  function  defined  in  Eq  (5)  for  pres¬ 
surization  loading  of  the  grain,  SIGPR 
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and  all  other  nomenclature  ia  the  same  aa  before. 


<V 


+ 1 


[<l>  ■*][  1-‘'c]<b»V 


(l+y)(h)(Ec) 


(5) 


A  complete  development  of  these  relatione  is  given  in  Reference  17. 


When  calculating  strains  in  the  CP  portion  of  a  grain,  a  =  R2  = 
radius  of  CP  port,  and  imposed  hoop  strains  will  be  found  directly  from 
Eqs  (3)  and  (4)  for  thermal  (EPT)  and  pressurization  (EPP)  conditions, 
respectively. 


When  calculating  strains  in  the  valley  of  star  or  slotted  portions 
of  the  grain,  a  =  radius  to  "bottom"  of  valley  for  use  in  Eqs  (3)  and  (4). 

Then  imposed  strains  are  found  by 

EPTS  =  Thermal  strain  in  valley  =  (K)  (EPT)  (6) 

EPPS  =  Pressurization  strain  in  valley  =  (K)(EPP)  (7) 


where  K  =  concentration  factor.  Figure  77  shows  a  section  of  a  typical  star 
geometry.  The  general  form  of  the  equation  for  the  concentration  factor  is 


K  = 


H 


(8) 


in  which  X  =  b/a  and  H  depends  on  the  geometry  of  the  star.  When  the 
angle  /9is  zero,  H  is  given  by 


H  =  N 


-1/3 


(9) 


where 


N 

P 

P 


Number  of  star  points  (2£  Ns  8),  NSLOTS 

Fillet  radius  (in)  between  star  point  and  web,  R4 

Included  angle  (deg)  of  valley  in  star  or  finocyl 
grain,  BETA 


There  are  three  basic  types  of  "star"  configurations  of  interest 
here.  They  are  the  finocyl,  the  star  and  the  wagon  wheel,  shown  in 
Figures  78,  79  and  80,  respectively.  They  will  now  be  discussed  individ¬ 
ually. 
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Finocyl 


The  angle  (3  of  the  slotted  tube  section  in  Figure  78  will  always  be 
less  than  40*  because  of  design  practice.  Figure  81  shows  the  dependance 
of  H  on  p  „  For  p  less  than  40°  the  effect  on  H  is  small,  so  for  this  con¬ 
figuration  Eq  (9)  will  be  used  to  calculate  H. 

As  an  example,  consider  a  finocyl  geometry  (Figure  78)  with  a  =  4, 
b  =  5,  d  =  2p  =  1,  N  =  5  and  p  =  21  *.  In  this  case  H  is  given  by  Eq  (7), 

H  =  1 1.  6789.  Then  X  =  1.  25  and  K  =  2,  10  from  Eq  (8). 

Star 


In  the  star  configuration  of  Figure  79  the  angle  8  will  not  necessarily 
be  less  than  40°.  Furthermore,  the  width  of  the  end  of  the  star  valley,  d, 
may  be  larger  than  twice  the  corner  radius,  p.  Modifications  to  Eq  (9) 
are  in  order  for  either  or  both  situations.  Consider  first  the  case  in  which 
P  is  larger  than  40°  and  d  =  2p.  In  th?s  case  H  is  given  by 


H 


H* 


L 


1 ±i_V5L 

i  +  2  -yiY 


(10) 


in  which  H*  is  selected  from  the  graphs  in  Figure  81  for  the  proper  g  and 
a/b  ratio. 

In  the  case  of  d>  2p  and  gs  40°,  the  H*  is  selected  from  Figure  82 
which  gives  the  dependance  of  H  on  the  d/2p  ratio,  and  Eq  (10)  is  used  to 
calculate  H. 

When  P  >  40°  and  d/2p  >  1  concurrently,  H*  must  reflect  both  facts. 
Define  a  factor,  F,  as 

F  =  H*  (d/2p)/H*(l)  (11) 

where  H*  (d/2p)  =  Photoelastic  parameter  obtained  from  Figure  82  at  the 
current  value  of  d/2p,  HSTAR2 

H*(l)  =  Photoelastic  parameter  obtained  from  Figure  82  at 

d/2p  =  1,  HSTARl 

Then  H  is  calculated  using 


H 


N 


-2/3 


=  FH*  (-£) 


1+2  Va  /  p 
i  +  2  -yi2 


2 


(12) 
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where  H*  =  photoelastic  stress  parameter,  HSTAR,and  is  obtained  from 
Figure  81  using  the  curr  s.pt  value  of  0. 

Consider  another  example  with  d  =  1.6,  p  =  0.  5  and  the  other 
dimensions  the  same  as  in  the  previous  example.  H  is  calculated  using 
Eq  (10)  with  H*  obtained  from  Figure  82.  H*  =  14.2,  H  =  11.  0683  and 
K  =  1.99.  Notice  that  widening  the  star  width,  d,  lowers  the  concentration 
factor,  K. 

W agon  Wheel 

The  effect  of  the  wagon  wheel  (Figure  80)  included  angle  a  on  II  is 
shown  to  be  pronounced  in  Figure  83.  Furthermore,  d/2p  will  generally 
be  greater  than  one,  so  the  factor  F  in  Eq  (11)  must  be  utilized  with  H* 
given  in  Figure  83  for  the  a  of  the  slot  to  determine  H  from 


H 


FH*  (^) 


2/3 


1  + 

[i*2  vsl 

i  +  2-yiz 

M. 

i  +  2  ys 

(13) 


Comparison  vith  Propellant  Capability 


Nominal  strain  endurance  is  furnished  by  the  user  as  a  fixed  value 
or  can  be  calculated  with  a  user- supplied  model  (see  User  Models  section 
of  this  volume).  Design  strain  endurance  is  derived  from  the  nominal 
value  by  accounting  for  statistical  variations  in  the  nominal  and  the  degra¬ 
dation  due  to  aging 


whei  e 


SEDES 

SEDES 

SENOM 

CVPS 

AGE 


SENOM  [  1  -  (3)(CVPS)  ]  [  1  -  AGE]  (14) 

Design  strain  endurance  (in/in) 

Nominal  strain  endurance  (in/in) 

Coefficient  of  variation  of  strain  endurance 
(%  x  0.01) 

Fraction  of  propellant  strain  endurance  las  t 
as  a  result  of  aging. 


After  applying  a  factor  of  safety  to  the  predicted  strain,  appropriate 
margins  of  safety  are  calculated 


MSP 

MS  PS 


SEDES 

(FSPSKEFT) 

SEDES 

(FSPS)(EPTS) 


(15) 

(16) 
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Photoelastic  Stress  Parameter,  H* 


From;  CPIA  Publication  214,  pg.  3.  66 


Positive  Wedge  Angle,  a  (deg) 


Figure  83. 


Positive  Wedge  Angle  Test  Results,  N  = 


i 


,  a/p  =  8 


6 


where 

MSP 

Margin  of  safety  for  storage  thermal  strain 
in  CP  portion  of  grain 

MBPS 

Margin  of  safety  for  storage  thermal  strain  in 
slotted  portion  of  grain 

FSPS 

Factor  of  safety  for  propellant  thermal  strain 

If  the  margins  of  safety  are  less  than  zero,  penalties  are  calculated  and 
summed  into  the  overall  penalty  function 

If  MSP<  0.0,  OBJPS  =  (MSP)2  105  (17) 

If  MS  PS  <  0.  0,  OBJPSS  =  (MS  PS)2  1  05  (1  8) 

Strains  due  to  ignition  pressurization  are  compared  with  a  user  -  supplied 
maximum  limit  (EPPMAX)  and  appropriate  penalties  calculated 

If  EPP  >  EPPMAX,  OBEPP  =  (EPP-EPPMAX)2  1  05  (19) 

If  EPPS  >  EPPMAX,  OBEPPS  =  (EPPS- EPPMAX)2  105  (Z0) 

Algorithm  Summary 

The  following  outlines  give  the  algorithms  used  to  evaluate  propel¬ 
lant  structural  integrity,  referring  to  the  equations  numbered  as  above  ancl 
the  accompanying  figures,  and  employing  the  nomenclature  used  in  the 
individual  grain  setup  subroutines  (see  Volume  II  for  definitions).  Also 
shown  is  the  internal  code  nomenclature  from  the  ballistic  simulation  module. 


Grain  Type  1  (Star) 

RZ  =  RFA1-TAUW1  =  PLANE  (1,4)  -  PLANE  (1,6) 

R4  =  R5A1  =  PLANE  (1,10) 

RF  =  RFA1  =  PLANE  (1,4) 

DOVERP  =  LITD/ (Z)(R4),  where  LITD  calculated  in  SETUPl 
subroutine 


If  BETA  s  40.  and  DOVERP  = 


SIGTHM 

SIGPR 

EPT 

EPP 

K 

EPTS 

EPPS 


from 

from 

from 

from 

from 

from 

from 


Eq  (Z) 
Eq  (5) 
Eq  (3) 
Eq  (4) 
Eq  (8) 
Eq  (6) 
Eq  (7) 
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If  BETA  s  40  and  DOVERP  >1.0 


HSTAR 

from  Figure  82  at  DOVERP  and  R2/RF 

SIGTHM 

from  Eq  (2) 

SIGPR 

from  Eq  (5) 

EPT 

from  Eq  (3) 

EPP 

from  Eq  (4) 

H 

from  Eq  (10) 

K 

from  Eq  (8) 

EPTS 

from  Eq  (6) 

EPPS 

from  Eq  (7) 

If  BETA  >40  and  DOVERP  =  1.0 

HSTAR 

from  Figure  81  at  BETA  and  R2/RF 

SIGTHM 

from  Eq  (2) 

SiGPR 

from  Eq  (5) 

EPT 

from  Eq  (3) 

EPP 

from  Eq  (4) 

H 

from  Eq  (10) 

K 

from  Eq  (8) 

EPTS 

from  Eq  (6) 

EPPS 

from  Eq  (7) 

If  BETA  >40  and  DOVERP  >1.0 

HSTAR2 

from  Figure  82  at  DOVERP  and  R2/RF 

HSTARl 

from  Figure  82  at  DOVERP  -  1  and  R2/RF 

F 

from  Eq  (11) 

HSTAR 

from  Figure  81  at  BETA  and  R2/RF 

H 

from  Eq  (12) 

SIGTHM 

from  Eq  (2) 

SIGPR 

from  Eq  (5) 

EPT 

from  Eq  (3) 

EPP 

from  Eq  (4) 

K 

from  Eq  (8) 

EPTS 

from  Eq  (6) 

EPPS 

from  Eq  (7) 

Grain  Type  2  (Wagon  Wheel) 

R2  =  RFA1-TAUW1  =  PLANE  (1,4)  -  PLANE  (1,6) 

R4  =  R5A1  =  PLANE  (1,10) 

RF  =  RFA1  =  PLANE  (1,4) 

DOVERP  =  LITD/(2)(R4),  where  LITD  calculated  in  SETUP2 
subroutine 

HSTAR2  from  Figure  82  at  DOVERP  and  R2/RF 
HSTAR  1  from  Figure  82  at  DOVERP  =  1  and  R2/RF 
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F  from  Eq  (11)  ,  _.,/DTP 

HSTAR  from  Figure  83  at  ALPHA  and  R2/RF 

H  from  Eq  (13) 

SIGTHM  from  Eq  (2) 

SIGPR  from  Eq  (5) 

EPT  from  Eq  (3) 

EPP  from  Eq  (4) 

K  from  Eq  (8) 

EPTS  from  Eq  (6) 

EPPS  from  Eq  (7) 

Grain  Type  3  (Finocyl) 

For  Slotted  Region 


R2  =  R5A1  =  PLANE  (4,9) 
RF  =  RFA1  =  PLANE  (1,  9) 
R4  =  R4A1  =  PLANE  (1,8) 
SIGTHM  from  Eq  (2) 

SIGPR  from  Eq  (5) 

EPT  from  Eq  (3) 

EPP  from  Eq  (4) 

H  from  Eq  (9) 

K  from  Eq  (8) 

EPTS  from  Eq  (6) 

EPPS  from  Eq  (7) 


For  CP  Region 

R2  =  R2A10  =  PLANE  (10,6) 
RF  =  RFA3  =  PLANE  (3,4) 
SIGTHM  from  Eq  (2) 

SIGPR  from  Eq  (5) 

EPT  from  Eq  (3) 

EPP  from  Eq  (4) 

Grain  Type  4  (Conocy;j 

R2  =  R2A3  =  PLANE  (3,  6) 
RF  =  RFA3  =  PLANE  (3,  4) 
SIGTHM  from  Eq  (2) 

SIGPR  from  Eq  (5) 

EPP  from  Eq  (4) 


Grain  Type  5  (CP) 

R2  =  R2A5  =  PLANE  (5,6) 

RF  =  RFA5  =  PLANE  (5,4) 

SIGTHM  from  Eq  (2) 

SIGPR  from  Eq  (5) 

EPT  from  Eq  (3) 

EPP  from  Eq  (4) 

Figures  81,  82  and  83  are  described  in  the  code  as  Tables,  each  as  an 
individual  subroutine;  these  subroutines  are  named  FIGS,  FIG6,  FIG7, 
respectively,  to  correspond  to  a  figure  numbering  system  of  the  Reference 
25  report. 

A  two-dimensional  plane-strain  model  is  used  to  calculate  propellant 
strain  due  to  low-temperature  storage  and  ignition  pressurization.  Such  a 
model  accurately  describes  the  propellant  behavior  at  a  point  mid-way  along 
the  grain  length  when  the  grain  length-to- diameter  ratio  (L/D)  is  equal  to 
or  greater  than  about  seven.  For  L/D  <  7,  or  for  locations  near  the  grain 
terminations,  the  plane-strain  models  give  very  conservative  predictions 
because  the  end  effects  (three-dimensional)  that  relieve  the  strain  are  not 
accounted  for  in  SPOC.  Strain  predicted  for  a  propellant  valley  or  slot  will 
also  be  conservative  near  the  ends  or  for  short  slots. 

The  propellant  structural  analysis  is  not  conservative  at  the  hinge 
points  of  stress  relief  flaps  and  at  the  transition  between  propellant  slots 
and  CP  regions.  Both  of  these  areas  represent  highly  three-dimensional 
conditions  that  are  not  amenable  to  preliminary  design  calculations  used  in 
SPOC.  Consequently,  there  is  the  inherent  assumption  that  the  bore 
conditions  are  the  critical  locations.  Provisions  have  been  made  to  include 
volume  and  weight  allowances  for  stress  relief  boots  ellipsoidal  closures, 
even  though  their  final  configuration  is  dependent  on  more  detailed  analyses. 
The  transition  section  between  slots  and  cylindrical  port  may  require  a 
special  configuration  to  limit  imposed  strains;  another  way  to  achieve  the 
same  results  is  to  specify  about  7  degrees  as  the  angle  on  the  side  of  the 
slot  (ALPHA1)  of  a  finocyl  grain  (Type  3). 

Thermal  strain  in  the  propellant  due  to  low-temperature  storage 
is  compared  with  design  strain  endurance  (nominal  strain  endurance 
red  uced  for  mix-to-mix  Variations  and  aging  degradation).  Strain  induced 
by  ignition  pressurization  is  compared  with  a  user-input  maximum  limit. 

This  latter  limit  should  be  derived  from  tests  that  measure  strain  capability 
at  rapid  strain  rate  (to  simulate  ignidon  pressurization  on  test  specimens 
conditional  to  the  design  low-temperature  and  already  strained  to  the  level 
that  will  be  induced  by  low  temperature  storage. 
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TRAJECTORY  SIMULATION 


Analytical  Relationships 

The  trajectory  simulation  is  based  upon  a  mathematical  model  of  the 
flight  dynamics  of  a  point-mass  missile  flying  a  two-dimensional  path  in  the 
altitude /range  plane  over  a  flat  earth.  Forces  modeled  are  restricted  to 
thrust,  drag,  and  weight.  The  ballistic  trajectory  restriction  assumes 
missile  orientation  such  that  lift  is  always  zero.  A  symmetrical  missile  is 
assumed,  resulting  in  angle  of  attack  always  being  zero.  At  zero  angle  of 
attack,  drag  and  axial  aerodynamic  force  are  equal.  There  is,  therefore, 
no  need  to  differentiate  between  the  two  common  force-accounting  systems 
(body  oriented  or  flight- path  oriented). 

The  time-dependency  of  thrust  and  propellant  weight  is  included  via 
the  output  of  the  motor  ballistic  simulation  subroutine.  Instantaneous  missile 
weight  is  taken  as  launch  weight  less  the  integral  of  motor  weight  discharge 
rate.  Variation  of  drag  with  Mach  number  is  described  through  a  user¬ 
generated  input  table.  Provision  is  made  for  separate  aerodynamic  data  to 
be  input  for  power-on  and  power-off  phases.  The  variation  in  atmospheric 
properties  with  altitude  is  modeled  from  the  1959  ARDC  STD  Atmosphere  and 
the  MIL-STD-210A  Tropical,  Polar,  Hot,  and  Cold  Atmospheres.  The  user 
shall  choose  one  from  among  these  choices.  Required  integration  of  time- 
dependent  parameters  is  accomplished  using  a  fourth-order  Runge-Kutta 
procedure. 

It  is  assumed  that  the  missile  will  be  air-launched  and  that  no  on- 
launcher  kinematics  will  be  included.  Launch  conditions  are  specified  by 
altitude,  Mach  number,  and  flight- path  angle.  A  pre-boost  glide  phase  is 
not  included,  i.  e.  ,  boost  ignition  occurs  at  the  instant  of  launch.  Provision 
is  included  for  two  flight  phases;  (1)  rocket  thrusting  and  (2)  post-burnout 
glide.  The  user  has  the  capability  of  terminating  the  trajectory  pimulation 
by  his  command.  .  He  may  specify  a  termination  upon  achieving  a  selected 
value  for  (1)  time  of  flight  from  launch,  (Z)  time  of  flight  after  boost  burnout, 
(3)  slant  range,  (4)  horizontal  range,  (5)  altitude  (either  approaching  from 
above  or  below),  (6)  missile  Mach  number,  (7)  missile  velocity,  (8)  flight  path 
angle,  (9)  missile^acceleration  along  the  flight  path,  or  (10)  range  along  the 
flight  path.  Also  termir  ±tion  may  be  commanded  upon  ground  impact  or 
boost  burnout.  The  boost  phase  is  always  completed  unless  ground  impact 
occurs  first.  The  termination  commands  apply  only  at  or  after  boost  burnout. 
Termination  upon  ground  impact  will  be  automatic  should  this  occur  before 
user- commanded  termination. 

The  trajectory  analytics  are  shown  on  the  following  pages. 
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NOMENCLATURE 


Aerodynamic  reference  area  for  the  miBBile  in  ft* 
Ambient  sonic  velocity  in  ft/sec. 

Missile  drag  in  pounds. 

Gravitational  constant  in  ft/sec^. 

Missile  flight  Mach  number. 

2 

Ambient  pressure  in  pounds/in  . 

Freestream  dynamic  pressure  in  Ibs/ft*" 

Time  in  seconds. 

Thrust  in  pounds. 

Missile  velocity  in  ft/sec. 

Missile  weight  in  pounds. 

Range  in  feet. 

Altitude  in  feet. 

Ambient  air  density  in  —  • 


(7) 

(8) 
(9) 


Flight  path  angle  in  degrees. 
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.  z 

""■ertical  acceleration  in  ft/aec  . 

2 

Horizontal  acceleration  in  ft/aec  . 

Nozzle  exit  area  in  in^. 

Mi8sile  drag  coefficient. 

Vacuum  thrust  in  pounds. 

Missile  velocity  at  launch  in  ft/sec. 

Horizontal  component  of  missile  velocity  in  ft/sea. 

Vertical  component  of  missile  velocity  in  ft/sec. 

Missile  weight  at  launch  in  pounds. 

Missile  range  at  launch  in  feet. 

Missile  altitude  at  launch  in  feet. 

Missile  flight  path  angle  at  launch  in  degrees. 

Propellant  weight  flowrate  in  pounds/sec. 

Slant  range  from  launch  in  feet. 

Horizontal  component  of  missile  velocity  at  launch  in 
ft/ sec. 

Vertical  component  of  missile  velocity  at  launch  in  ft/sec. 
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Execution  Logic 


If  the  user  select*  FTRAJ=T  in  namelist  CONTRL,  a  trajectory  simu¬ 
lation  will  be  performed.  Subroutine  TRAJIN  reads  the  user  inputs  and 
digests  tho  input  data  *lnd  subroutine  TRAJ  performs  the  simulation.  TRAJ 
is  called  by  TRAJIN.  If  the  user  is  analyzing  a  two-temperature  problem, 
and  his  inputs  require  running  the  trajectory  simulation  at  both  low  fcnd 
high  temperatures,  then  TRAJIN  will  be  called  only  on  the  low  tempera¬ 
ture  trajectory  simulation  (i.  e.  ,  using  the  low- temperature  thrust,  history  K 
The  variable  DELTAV  (change  in  velocity)  is  used  as  a  flag,  and  a  non-zer  > 
value  indicates  that  TRAJIN  (and  the  subsequent  low-temperature  trajectory 
simulation)  have  been  run.  Figure  84  shows  the  trajectory  decision  logic. 

If  any  of  the  requirements  associated  with  low-temperature  trajectory 
analysis  are  different  from  their  default  values  (DE  LVRQ  >  1,  TTTRQ  <.  .  99E6, 
or  VTRQ  >  1),  TRAJIN  will  be  called  with  the  low-temperature  thrust  history 
tables.  If  this  test  falls,  NOLO  will  be  set  to  1,  which  indicates  no  low- 
temperature  requirement  has  been  specified.  If  the  requirement  associated 
with  high- temperature  trajectory  analysis  is  different  from  its  default  value 
(ACLIM<  .  99E6),  and  the  ballistics  simulation  was  performed  at  only  one 
temperature,  the  flag  DELTAV  indicates  that  TRAJIN  should  be 
callbd.  If  ACLIM  <  .  99E6  and  ballistics  were  run  at  both  low  and  high  temper¬ 
atures,  then  the  high-temperature  thrust  tables  are  used  to  run  the  trajectory 
simulation.  DELTAV  indicates  whether  TRAJIN  or  TRAJ  should  be  called. 

If  the  test  fails,  then  NOHI  is  set  equal  to  1,  which  indicates  no  high  temper¬ 
ature  requirement. 

If  the  user  inputs  no  low  or  no  high  temperature  requirements,  a 
trajectory  simulation  will  be  performed  only  on  the  first  and  last  passes 
through  the  program  if  FTRAJ=T  is  selected  by  the  user. 

Ideal  drag-free  burnout  velocity  and  axial  acceleration  are  calculated 
in  subroutine  FLT  if  FTRAJ=F  (its  default  value).  Velocity  at  launch  is 
assumed  equal  to  zero.  Burnout  velocity  is  calculated  at  both  temperatures 
in  a  two- temperature  problem,  but  the  low  temperature  value  is  compared 
with  the  requirement  (DELVRQ).  Axial  acceleration  is  calculated  at  high- 
temperature  for  comparison  with  its  requirement  (ACLIM). 
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Integration  Technique 


For  those  instances  in  which  it  will  be  necessary  to  integrate  time- 
dependent  variables,  a  fourth-order  Runge-Kutta*  numerical  integration 
procedure  will  be  employed.  This  method  is  widely  used  because  of  its 
many  desirable  characteristics,  including:  (1)  it  is  a  single-Btep  method 
and  is  therefore  self-starting,  (2)  its  accuracy  is  great  relative  to  the 
independent  variable  step-size,  (3)  independent  variable  step-size  may  be 
changed  at  any  time  without  affecting  previous  computations,  and  (4)  ihe 
method  is  readily  adapted  to  systems  of  simultaneous  equations  where  integra¬ 
tion  in  parallel  is  required. 

The  necessity  for  numerical  integration  in  the  SPOC  computer  code 
arises  as  a  result  of  the  need  to  perform  trajectory  simulations.  The  rocket 
motor  weight:  discharge  rate  must  be  integrated  to  determine  instantaneous 
missile  weight.  The  missile  acceleration  must  be  integrated  once  to  deter¬ 
mine  velocity  and  again  to  determine  displacement.  The  nature  of  these 
integrations  differ  substantially  and  therefore  the  application  of  the  Runge- 
Kutta  integration  procedure  to  each  will  be  treated  seoaratelv. 

Rocket  Motor  Weight  Discharge  Rate 

The  rocket  motor  internal  ballistics  subroutine  will  produce  a  schedule 
of  rocket  motor  weight  discharge  rate  as  a  function  of  time.  This  complete 
schedule  of  motor  weight  discharge  rate  with  time  is  known  prior  to  any 
trajectory  integration  being  performed.  The  form  of  these  data  is  a  simple 
stored  table  of  motor  discharge  rate  versus  time.  The  fact  that  the  rocket 
motor  weight  discharge  rate  is  independent  of  missile  flight  dynamics 
allows  these  data  to  be  integrated  separately  from  (i.  e.  ,  not  parallel  to) 
missile  flight  dynamics  and  therefore  special  simplifications  apply  to  the 
general  Runge- Kutta  numerical  integration  procedure.  In  this  special 
case  where  the  derivative  of  rocket  motor  weight  discharged  with  respect 
to  time  is  a  function  of  time  only,  the  general  fourth- order  Runge- Kutta 
numerical  integration  procedure  defaults  to  the  same  scheme  known  as 
Simpson's  Rule.  For  this  special  case, 


dWd 

—  =  *d  =  f<t)  (10) 

where  wd  is  the  time-derivative  of  rocket  rro  tor  weight  discharged  and  is  the 
function  to  be  integrated.  The  result  of  this  integration  is  the  decrease  in 
rocket  motor  weight  (also  the  decrease  in  missile  weight)  to  a  given  point  in 


*A  Basic  Course  In  Numerial  Methods;  Ralph  E.  Ekstrom;  reprinted  from 
Machine  Design;  October  2  6,  1967  through  June  20,  1968;  Penton  Publishing 
Co.;  Cleveland,  Ohio  44113 
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time.  In  stepping  through  the  numerical  integration  process  from  time  t 
to  time  =  t  +  At, 


w  =  w  +  (At)  (w  )  (11) 

t  +  £t  t 

where 

w  is  the  rocket  motor  weight  discharged  to  time  -  t  +  At 

dt  + 

w,  is  the  rocket  motor  weight  discharged  to  time  -  t 
d 


At  is  the  numerical  integration  step  size 

w  is  the  average  value  of  w  between  time  =  t  and  time  =  t  +  At 


The  numerical  integration  procedure  is  concerned  with  the  evaluation  of  the 
term  (Ath^d)  in  equation  (ll).  For  this  special  case,  the  fourth-order 
Runge-Kutta  method  estimates  this  term  to  be 
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where 


w  is  the  rocket  motor  weight  discharge  rate  at  time  -  t 
dt 

w  is  the  rocket  motor  weight  discharge  rate  at  time  =  t  +  ^ 

dt  +  ^ 

+  2 


w,  is  the  rocket  motor  weight  discharge  rate  at  time  =  t  +  A* 

d 

t  +  At 
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The  mechanics  of  the  numerical  integration  procedure  consists  of  a 
cyclic  repetition  of  the  fcllowing  steps: 


(1)  The  results  from  the  rocket  motor  internal  ballistics  subroutine 
make  known  the  value  of  w<j  at  time  =  t.  Select  a  value 
for  At  And  compute  kg  from  equation  (13). 


(2)  The  results  from  the  rocket  motor  internal  ballistics 

subroutine  makes  known  the  value  of  w^  at  time  =  t  t  &t. 

Then  compute  k^  from  equation  (14).  2 


(3)  The  results  from  the  rocket  motor  internal  ballistics  subroutine 
makes  known  the  value  of  W(jt  at  time  =  t  +  A**  Then 

compute  k^  from  equation  (15). 


(4) 


Knowing  kg,  kj,  and  k?, 
and  (12), 


solve  for 


from  equations  (11) 


(5)  Reset  time  t  to  t  +  At  and  return  to  Step  (1).  Continue  cycling 
until  time  becomes  equal  to  motor  burn  time  determined  from 
the  rocket  motor  internal  ballistics  subroutine. 


The  above  procedure  defines  rocket  motor  weight  discharged  as  a 
function  of  time  from  ignition  to  burnout. 


Missile  Acceleration 

The  rocket  motor  infernal  ballistics  subroutine  will  produce  a  schedule 
of  vacuum  thrust  as  a  function  of  time.  This  complete  history  of  vacuum 
thrust  is  known  prior  to  any  miesile  flight  dynamics  integration  oeing  per¬ 
formed.  The  complete  history  of  motor  weight  discharged  during  the  burn 
is  also  known,  this  information  having  been  determined  rs  a  result  of  the 
previously  described  integration.  It  follow*  that  the  missile  weight  is  known 
a:  any  instant  of  time.  The  integration  of  missile  acceleration  is  then  per¬ 
formed  in  the  following  manner; 

.  f  (t.  V  ,  V  )  (16) 

*  dt  7 


where 


a  is  the  missile  horizontal  acceleration 
x 

a  is  the  missile  vertical  acceleration 

y 

Vx  is  the  missile  horizontal  velocity 
Vy  is  the  missile  vertical  velocity 


dV 

_ }_  - 

dt 


j  («,  V. 


(17) 
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f  (  )  ia  the  function  defining  missile  horizontal  acceleration  and  ia  the 
system  of  equations  numbered  (2),  (3),  (4),  (5),  (6),  (7),  (8), 
and  (9)  in  the  TRAJECTORY  section. 

j  (  )  is  the  function  defining  miaaile  vertical  acceleration  and  is  the  system 
of  equations  numbered  (1),  (3),  (4),  (5),  (6),  (7),  (81,  and  (9) 
in  the  TRAJECTORY  section. 

t  is  the  independent  variable  time 

In  stepping  through  the  numerical  integration  process  from  time  =  t 
to  time  =  t  +  At. 


V 

*t  +  At 


Vx  +  (At)(a  ) 

xt  x 
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V  =  V 

Vt  +  At  yt 


+  (At) (a  ) 

y 


(19) 


where 


V 

*t  +  £t 


is  missile  horizontal  velocity  at  time  stt^t 


V 

Yt  +  £t 


is  missile  vertical  velocity  at  time  =  t  +  At 


V„  is  missile  horizontal  velocity  at  time  =  t 
*t 


Vyt  is  missile  vertical  velocity  at  time  =  t 

At  is  the  numerical  Integration  step  size 

is  the  average  of  a^  between  time  =  t  and  time  -  t  + 

The  numerical  integration  procedure  is  concerned  with  the  evaluation 
of  the  terms  (AOf*x)  Hnt*  (At)(»y)  in  equations  (18)  and  (19)  respectively.  The 
fourth-order  Runge-Kotta  numerical  integration  procedure  estimates  these 
terms  to  be 

(At)(*x)  =  (1/M(k0  +  2kj  +  2k2  x  k3j  (20)  (At>(*y>  =  (1/6)(L0  +  2Lj  +  2L,  +  L3) 

(21) 

k0  -  (At)  «  (t.  VXfc,  Vyt)  (22)  Lq  =  (At)  j  (t.  VXt.  Vy)  (23) 


190 


kl  =  (tt)U t+f  .  VX(  +  Vy[  + 


k  L 

Lj  =  (At)j(t+AS.  vx^  +  ~  ,  v  +  —■) 


At  ^1  ^*1 

*2  =  (At)f(t+^  vXt  i  T.  V  +  Y) 


L2  -  v,t  +  j;  vy  +  ii| 


k3  =  (At)f(t  +  At.  vx  +k2,  V  +  l2) 


L3  -  (At)  j  (t  +  At,  VXfc  +  k2,  Vy^  +  L2) 


where 

f  (t,  Vx  ,  Vy  )  is  the  function  f  (  )  evaluated  at  time  =  t,  missile 

1  t  horizontal  velocity  =  V  ,  and  missile  vertical  veloc  ity 


j  (t,  Vx  ,  Vy  )  is  the  function  j  (  )  evaluated  at  time  =  t,  missile 

t  t  horizontal  velocity  =  V  ,  and  missile  verticle  velocity  - 


At  *0  L0 

f  (t  +  “r- ,  V  t  —  ,  V  +  ~r~)  is  the  function  f  (  )  evaluated  at  time 
£  t  4  £- 

t  ^  A",  missile  horizontal  velocity  =  V  +  —  , 

Lo 

and  missile  vertical  velocity  =  V  +  -r — 

yt  2 
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k  L. 

j  (t  +  AL.  t  y  +  ,  V  +  -y )  is  the  function  j  (  )  evaluated  at  time  = 

2  xt  2  yt  t  +  Ai.  mieelle  horizontal  velocity  - 

2 

k. 

V  +  ~r~  ,  and  missile  vertical  velocity  = 
xt  2 


Vyt  +  T 


f  (t  +  ^  ,  V  +  — .  V  +  —  )  is  the  function  f  (  )  evaluated  at  time  =  t  +  •*- 
Xt  yt  missile  horiaontal  velocity  =  V  +  kl  ,  and 

*L  2 

missile  vertical  velocity  =  Vv  +  1 

k.  L,  . t 

j  (t  +  Ap-  ,  vx  +  —  ,  Vy  +  -~)  is  the  function  j  (  )  evaluated  at  time  t  + 

missile  horizontal  velocity  =  V  +  -r~,  and 

t  2 3 

Lj 

missile  vertical  velocity  =  Vv  +  -r- 

yt  2 

f  (t  +  ^t,  Vx  +  k  ,  Vy  +  L,)  is  the  function  f  (  )  evaluated  at  time  =  t  +  frt, 

4  4  missile  horizontal  velocity  =  Vx  +  k^,  and 

missile  vertical  velocity  =  Vy  +* 

j  (t  +  At.  V  +  K  ,  V  +  Lz)  i«  the  function  j  (  )  evaluated  at  time  = 
t  L  yt  t  +  A1*  missile  horizontal  velocity  = 

V  +  K,,  and  missile  vertical  velocity  = 

Xt 

V  +  L, 
y  2 


At 

t  + 


The  mechanics  of  the  numerical  integration  procedure  consis’s  of  a 
cyclic  repetition  of  the  following  steps; 

(1)  The  results  of  the  rocket  motor  internal  ballistics  subroutine 
determine  the  value  of  vacuum  thrust  at  any  point  in  time. 

The  previously  performed  integration  of  motor  weight  dis¬ 
charge  rate  determines  missile  weight  at  any  instant  in  time. 
Selection  of  a  value  for  A*  will  now  permit  the  determination 
of  kQ  and  LQ  from  equations  (22)  and  (23)  respectively. 

(2)  Knowing  kQ  and  LQ,  compute  kj  and  Lj  from  equations  (24) 
and  (25). 

(3)  Knowing  kj  and  Lj,  compute  k£  and  from  equations  (26) 
and  (27V 
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(4)  Knowing  k,  and  L?,  compute  k  and  L  from  equation*  (28) 

and  (29).  '  3  3 

(5)  Knowing  kQ,  k^  k^  k3>  LQ,  Lj,  L.,,  and  L 3  compute  (£t) 

(a  )  from  equation  (20)  and  (At)(a  )  from  equatior  (21).  Knowing 
x  y 

(At)(*x)  and  (At)(ay),  solve  for  VXt  +  and  Vyfc  +  ^  from 
equations  (18)  and  (19). 

(6)  Reset  t  to  t  +  £t  and  return  to  step  (1).  Continue  cycling  until 
a  trajectory  termination  command  is  encountered.  The  result 
will  be  a  history  of  horizontal  and  vertical  velocity  throughout 
the  missile  time  of  flight. 

Missile  Velocity 

Once  the  integration  of  missile  acceleration  has  been  performed, 
a  history  of  missile  velocity  versus  time  is  known  throughout  the  entire  flight. 
Thus  missile  velocity  is  a  function  of  time  only  and  its  integration  will  be  a 
special  case  similar  to  that  previously  described  for  motor  weight  discharge 
rate.  As  was  noted  previously,  this  special  case  of  the  fourth-order  Runge- 
Kutta  procedure  defaults  to  Simpson's  Rule.  The  integration  of  missile  veloc¬ 
ity  may  be  performed  either  (1)  in  series  with  the  integration  of  missile 
acceleration  (that  is  one  after  the  other)  or  (2)  in  parallel  with  the  integra¬ 
tion  of  missile  acceleration  (that  is  both  at  the  same  time).  SPOC  performs 
acceleration  and  velocity  integration  in  parallel.  The  integration  of  missile 
velocity  will  define  missile  translation  and  is  performed  in  the  following 
manner: 

3?  =  vx  =  h  <t)  po>  fj:  *  vy  *  b  <‘>  on 


where 

V  is  the  missile  horizontal  velocity 
x 

V  is  the  missile  vertical  velocity 

y 

h  (t)  is  the  furction  defining  missile  horizontal  velocity  as  a  function 
of  time  ?nd  is  the  result  of  the  previously  described  integration 
of  missile  acceleration, 

b  (t)  is  the  function  defining  missile  vertical  velocity  as  a  function  of 
time  and  is  the  result  of  the  previously  described  integration 
of  missile  acceleration. 

In  stepping  through  the  numerical  integration  process  from  time  =  t 

to  time  =  t  +  A  t. 
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X  =  X  +  (At)(V  )  (32) 

t  +  A*  t  “  *  '  ' 


t  +  A1 


Yt  +  (At)(Vy) 


(33) 


where 


X 

Y 

X 

Y 


t 

t 

t 

t 


+  is  missile  horizontal  range  at  time  =  t  +  At 

,  is  missile  altitude  at  time  =  t  +  At 
+  At 

is  missile  horizontal  range  at  time  =  t 
is  the  missile  altitude  at  time  =  t 


At  is  the  numerical  integration  step  size 

V  is  the  average  misaile  horizontal  velocity  between  time  =  t 
and  time  =  t  +  At 


V  is  the  average  missile  vertical  velocity  between  time  =  t  and 
^  time  =  t  +  At 

The  numerical  integration  process  is  concerned  with  the  evaluation 
of  the  terms  (At)(V  )  and  (At)(V  )  in  equations  (32)  and  (33)respectively.  For 
this  special  case,  &e  fourth-orier  Runge-Kutta  method  estimates  these  terms 
to  be 


(At)(Vx)  =  (1/6)(K0 

+  4Kj  +  K2) 

(34)  (At)(Vy)  =  (1  /6)(N0  +  4Nj  +  N2) 

Kq  =  (At)(VXt) 

(36) 

=  (At)(Vyt)  (37) 

K1  =  <*><%  + 

2 

(38) 

Nl  =  (AtHVyt+At)  (39) 

2 

K2  =  <*><%  +  At) 

(40) 

N2  *  <At)(Vy|  +  At)  (41) 

where 

Vx  is  the  missile  horizontal  velocity  at  time  -  t  (obtained  from  previous 
integration  of  acceleration) 

^Yt  is  the  missile  vertical  velocity  at  time  *  t  (obtained  from  previous 
integration  of  acceleration) 
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V  Li  the  missile  horizontal  velocity  at  time  =  t  +  ^  (obtained  from 

**  ^  ^  previous  integiation  of  acceleration) 

Vy  is  the  missile  vertical  velocity  at  time  =  t  +  (obtained  from 

previous  integration  of  acceleration) 


V*t  +  At 


is  the  missile  horizontal  velocity  at  time  =  t  +  At  (obtained  from 
previous  integration  of  acceleration) 


V 

yt  +  At 


is  the  missile  vertical  velocity  at  time  =  t  +  A*  (obtained  from 
previous  integration  of  acceleration) 


The  mechanics  of  the  numerical  integration  procedure  consists  of  a 
cyclic  repetition  of  the  following  steps: 

(1)  The  results  of  the  previous  integration  of  missile  acceleration 
determines  missile  horizontal  and  vertical  velocity  at  any  instant 
of  time.  Select  a  starting  time  and  &t.  Solve  for  and 

from  equations  (3b)  and  (37). 

(2)  Solve  lor  Kj  and  Nj  from  equations  (38)  and  (39). 

(3)  Solve  for  and  from  equations  (40)  and  (41). 

(4)  Knowing  KQ,  Kj,  K2,  NQ,  Nj.  and  N2>  solve  for  (AWVJ 
from  equation  (34)and  for  (^t){V  )  from  equation 

(5)  Solve  for  X,  ,  . .  from  equation  (32)  and  for  Y.  from 

t  +  t  + 

equation  (33). 

(b)  Reset  time  t  to  t  +  A*  and  return  to  step  (1).  Continue  cycling 
until  a  trajectory  termination  command  is  encountered.  The 
result  will  be  a  complete  history  of  altitude  and  horizontal 
range  throughout  the  missile  flight. 


Integration  Step  Size  Determination 

In  using  any  numerical  integration  procedure,  the  allowable  error 
at  the  end  of  each  step  determines  the  interval  length.  If  the  interval  is 
smaller  than  necessary,  the  number  of  computational  cycles  will  be  unneces¬ 
sarily  great  and  excessive  computer  run  time  will  result.  If  the  interval  is 
too  large,  computational  accuracy  will  suffer.  The  desired  compromise  is 
to  select  an  interval  sufficiently  large  to  just  avoid  exceeding  needed  accuracy. 
The  needed  accuracy  now  becomea  a  judgement  criterion  and  must  be  either 
user  stated  or  implied  within  the  code.  Probably  a  separate  accuracy  criteria 
would  be  needed  for  each  integrated  parameter.  The  substantial  increase 
in  complexity  required  of  variable  step  alee  Integration  logit  was  »•«•» 
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to  be  justified.  SPOdC{l^^  * 

logic  dwini  P°w«red  Oi^tHttih.  P>im  urnQUt  co„t  flight  i.  ai.o 

££  ;;r:w  v^l  to  a  u.er-input  multiple  of  the  boost  powered 
flight  trajectory  step  size. 
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COST 


Two  option*  are  available  to  eatimate  coat*.  One  uses  the  Tri- 
Servicea  Coat  Study  rcaulta;  the  other  employ*  a  model  aupplled  by  the  uaer 
(see  Uaer  Model  aection  of  thia  volume).  Either  option  la  initiated  by  setting 
FCOST  =  T  in  the  nameliat  CONTRL.  If  the  Tri-Services  model  ia  desired, 
nothing  more  ia  required.  If  a  user  model  ia  supplied.  CSTMDL  =  T  muat 
also  be  included  in  CONTRL. 

Tri-5ervicea  Coat  Model 


This  option  employs  the  general  coat  relations  for  ftteel-case  motors 
developed  in  the  Tri-Services  Cost  Study  (Reference  18).  F  irst  unit  pro¬ 
duction  coat  (BFUCST)  of  the  basic  motor  ia  found  from 

BFUCST  =  (4.493 KWMOTOR)0,  1306  (ISP)0,  711  (1-MFf  828  (1) 

where 

WMOTOR  =  Total  motor  weight  (lbm) 

ISP  *  Delivered  specific  impulse.  ?0*F  (lbf-aec/lbm) 

MF  =  Motor  mass  fraction 

If  a  single-temperature  problem  la  being  run,  ISP7Q  calculated  for  that 
simulation  becomes  ISP.  If  a  two- temperature  problem  is  being  run,  ISP 
is  estimated  by 

In  (BP)  =  In  (ISPHI)  -  £  j  In  (ISPHI/ISPLO)  (2) 

Then  BFUCST  ia  adjusted  for  the  components  auch  as  igniters,  blast 
tubes,  etc.,  using  the  factors  in  Table  6.  The  code  uaer  must  select  from 
Table  6  those  multiplicative  factors  which  apply  to  his  particular  problem; 
the  product  of  these  individual  values  are  input  as  MULFAC.  Two  compo¬ 
nents  are  additive  factors;  igniter  and  safe-and-arm  device;  an  appropriate 
sum  is  input  as  ADDFAC.  Thus,  the  motor  first  unit  production  cost  (FUFCST) 
is 


FUPCST  =  (MULFAC  )(BFUCST)  +  ADDFAC 


(3) 


The  expected  level  of  production  is  used  to  adjust  the  first  unit  cost 
to  the  average  unit  cost.  Two  adjustments  are  needed;  production  rate 
(PRATE)  and  production  quantity  (PQUAN).  Learning  curves  of  96%  for  pro¬ 
duction  r  ate  and  94%  for  production  quantity  are  employed  for  this  adjustment 


PRATEF  = 
PQUANF  = 


1.314  [  (PRATE)"0,  0589 
1.  098  [  (PQUAN)"0,0893 


(1  /PRATE)] 
(1  /PQUAN)] 


(4) 

(5) 
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TABLE  6 


PRODUCTION  COST  FACTORS 
Adjustment _ Factor 


c  C*se  Attachments 

Forward  attachment  to  missile  1.11 

Launch  lugs,  aft  of  pressure  vessel  1.  07 

Launch  lugs,  on  pressure  vessel,  integral  1. 13 

Launch  lugs,  cn  pressure  vessel,  strap- on  1.  07 

Fin  attachment,  dovetail,  untapered  1.20 

Fin  attachment,  dovetail,  tapsred  1.30 

Fin  attachment,  folding  1.25 

Fin  clips,  fixed  1.07 

o  Blast  tube  i .  09 

o  Canted  nozzle  1.05 

o  Grain 

Dual  thrust,  single  grain  1.04 

Composite  smoky  1.00 

Composite  reduced  smoke  0.  98 

Double  base  smok>  1.28 

Double  base  minimum  smoke  1.44 

Dual  thrust,  dual  grains  1.  12 

High  burn  rate,  greate-*  than  1.  5  in/ sec  1.12 

High  burn  rate,  greater  than  3.  0  in/sec  1.30 

Free  standing  grain,  internal  burning  0.  80 

Free  standing  grain,  internal/ external  0.80 

burning 

o  Inert  slivers  1.07 

o  External  insulation  1.05 

o  Thrust  vector  control 

Liquid  injection  TVC  1.30 

Flexible  nozzle  1.40 

Hot  gas  bleed  1.40 

Warm  injection  and  jet  interaction  1,30 

Jet  vanes  1,35 

Jet  tabs  1.37 


Table  6  (Continued) 
PRODUCTION  COST  FACTORS 


_ Adjuitment _  Factor 

o  Boost /Sustain 

2;1  through  5;1  1.04 

6:1  1.06 

7;1  1.08 

8:1  1.10 

9:1  1.14 

10:1  1.20 

o  Pulse  Mode 

One  One  pulse  (two  grains)  1.33 

Two  pulse  (three  grains)  1. 64 

Four  pulse  (five  grains)  1.80 

o  Thermal  cookoff  1.04 

o  Thrust  termination 

Propellant  extinguishment  1.  09 

Thrust  reversal  1.13 

o  Wire  harness  1.09 

o  R I  filter  1.04 

o  Igniter  (additive)  $336 

o  Safe /arm  (additive) 

Manual  $462 

Remote  $2217 
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Thus  the  average  unit  production  coat  la 

AUPCST=(PRATEF)(:>QUANF){FUPCST)  (6) 

and  the  total  production  coat  becomea 

TPCST  =  (AUPCST)(PQUAN)  (7) 

Coata  for  development  (DEVCST),  pr^-flight  readineaa  testing  (PFRTCS) 
and  qualification  (QUALCS)  are  eatlmated  by 

DEVCST  =  ((0.  00379)(ISP)(BFUCST)  +  261.)  1000  (8) 

QUALCS  *((0.  000736)(ISP)(BFUCST)  +  51.)  1000  (9) 

PFRTCS  » ((0.000736HISP)(BFUCST)  +  51)  1000  (1C) 

Eq  (9)  and  (10)  are  Invoked  only  if  apecified  by  che  uae  (QUAL  *  T  and  PFRT  *  T, 
respectively). 

The  total  project  coat  (COST)  ia 

COST  *  DEVCST  +  PFRTCS  +  QUALCS  +  TPCST  (11) 


This  parameter  can  be  uaed  merely  aa  other  data  by  which  the  uaer  evaluatea 
the  design,  or  by  selecting  ICHOSE  •  1,  COST  becomea  thp  parameter  to 
be  minimised  by  the  optimizer. 

User  Model 


As  described  in  the  Uaer  Model  section,  the  user  builds  his  own  sub¬ 
routine,  including  whatever  common  blocks  that  are  necesaary  to  provide 
his  model  with  the  necesaary  input  from  other  parts  of  the  code.  The  param 
eter  COST  must  be  supplied  for  optimization. 


COMBUSTION  STABU  ITY  ANALYSIS 


Recognition  that  propulsion  system  mission  failure  or  degradation 
can  result  from  the  effects  of  combustion  instability  has  led  to  increasing 
emphasis  on  combustion  stability  as  a  design  parameter.  Accordingly, 
motor  stability  was  deemed  a  necessary  optimization  parameter.  The  one¬ 
dimensional  longitudinal  Standardized  Stability  Prediction  (SSP)  (Ref.  19  ) 
was  selected  for  use  in  this  program  due  to  its  general  acceptance  in  the 
combustion  community  and  general  agreement  with  experience  (Ref.  20  ). 

The  stability  analysis  module  is  a  version  of  SSP  that  has  been  modified 
to  reduce  execution  time  and  to  include  various  desirable  options  and  fea¬ 
tures.  The  code  is  tailored  for  use  as  an  optimization  tool  and  to  add  two 
additional  combustion  response  models. 

Except  as  noted  in  the  following  discussion,  the  same  philosophy, 
theory  and  general  coding  logic  are  used  as  in  SSP  and  will  not  be  repeated 
in  this  report.  The  vast  majority  of  variable  names  are  the  same  as  in 
SSP,  so  the  coding  will  largely  be  familiar  to  persons  experienced  with  the 
coding  in  SSP. 

Figure  85  is  a  diagram  of  the  general  organization  of  the  stability 
analysis  code  block.  Entry  to  this  stability  analysis  is  accomplished  by 
calling  subroutine  E488M2.  The  majority  of  the  required  data  is  trans¬ 
ferred  by  way  of  Common.  User  input  has  been  minimized  by  the  internal 
data  transfer  and  by  selecting  the  options  considered  most  appropriate 
to  tactical  rocket  motors.  The  only  direct  user  input  namelist  is  STAUIN, 
in  which  the  user  may  specify  the  number  of  modes  to  be  analyzed  and  the 
combustion  response  model  to  be  used.  Default  values  are  provided,  so 
even  STABIN  inputs  are  not  required. 

Figure  86  shows  the  locations  of  the  twenty  "sections"  used  to  describe 
the  motor  cavity  for  stability  analyses,  and  their  relationships  to  the  four 
'  stations"  in  the  head-end  and  nozzle-end  portions  of  the  grain  and  to  the 
fourteen  "planes"  used  to  define  the  center  portion  of  the  grain.  For  closure 
Type  1:  (1)  Section  ;  is  the  closure  itself  prior  to  ignition;  (2)  Section  2  is 
essentially  non-e?  ;  stant  at  ignition,  but  grows  in  length  as  propellant  is 
consumed  (Figure  86  shows  it  at  an  intermediate  position);  (3)  End  burning 
surface  at  the  aft  end  is  part  of  Section  18.  For  closure  Type  2:  (1)  Section 
1  length  does  not  char  ge  during  burning;  (2)  End-burning  surfaces  are  part 
of  Section  2  and  18  at  the  forward  and  aft  ends,  respectively;  (3)  Sections  2 
and  18  increase  in  length  and  Sections  3  and  17  decrease  in  length  as  pro¬ 
pellant  is  consumed  if  the  grain  ends  are  not  inhibited  (otherwise  their  length- 
stays  constant).  Table  7  gives  the  sources  of  all  data  used  in  the  stability 
analysis. 


201 


Catling 

Subroutine 


OURCE  OF  STABILITY  ANALYSIS  DATA 


2  ^  w.  »r  *r»  vri  *•••*•••'♦  un  •  ‘  o>  »n  i/>  kTi 

»M  ^•w.^^r-w'^.wOOOOOOOOOO—  oo-—  w  ^ 


irNr-^'£ir^00®®®00®00tAD0^irr'f^ 

^  ^  *r»  to  -*••*•••'■••*«%•*  »r»  »r»  »r»  »r» 

<—  —  wO  oo  o  o  o  o  a  ©  ©  —  ©o 


T  u 

IS  N  «  n  g,  ^  ^  <0  .  .  O  O'  .  .  .  ,  *  ,®nwNm»  « 

j  2  ^<i;CC.tOO**OOOOOONm<n«tt!Di 


««w?5  wCSw-. 

Z  Z  Z  Z  £  £  °.  °.  °.  ®  Z  Z  ZZ£- 

2555i^oooa3535-: 

0-  0-  Q,  Cl,  <  <  d«  Q,  (X  0* 


OOlNJN^OOiO^  -D  « 
N  N  N  N  M  m  ♦  t  1.1 


QQflQ 

.  _  W  W  W  W - - 


uuuu 

z  :i  z  z  — - - 


''CfiElSiSodoo  J  J  JJttttiSttH 
<  <  <  <  0.0.CL9. 


Q  Q  U  U 

W  W  5)  tn 

J  ^  J  ^  j—  j-  ^ 

"*  X~~A,‘o0.^oooo0Qooeo<°^o^T£ 
ait<"<t.odcio  X  X  o'  o'  o'  o'  —  o  -  o'  -  i  i 


U,  Q  Q 

_.  a.  w  u 


o  <  x  —  I  — . 

~l  i.^n£^a.^00000 

I  tl.  Q  j  <  C.  <  ”  o'  o'  o'  o'  <5 


fa. 

-  a, 

0000000000^.000*03 

OOOOOOOOOO^toOCwOi-4 


!*4  no  *n  -O  «*>  m  try  r**.  -\)  Tl  N  N  rO 


0 

Q  , 

•J 

Ui 

P4 

u 

U 

Q 

X 

u  0  J 

M  *J 

a.  x 

w 

Si 

CQ 

*0 

(Q 

W 

a 

<n 

a 

3 

J 

<  « 

<  I  J 
2  ce  x 

H  H 

H 

H 

H 

H 

H 

H 

H 

H  H 

H(-t- 

Table  7  (Contd. ) 


Table  7  (contd.  ) 

Source  of  Stability  Analysis  Data 


Table  7  (contd.  ) 

Sourc  e  of  Stability  Analysis  Data 


Table  7  (contd.  ) 

Source  of  Stability  Analysia  Data 
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PLANE(I,  79)  Combustion  gas  density  (lbm/cu  in)  at  Plane  I.  1  -  1.14 


Table  1  (contd.  ) 

Source  of  Stability  Analysis  Data 


Table  7  (contri,  ; 


E488M2 


This  is  the  control  subroutine.  It  wee  derived  from  SSP/MAIN. 
SSP  capabilities  and  subroutines  deleted  (because  not  needed)  were 


SSP/TABUL 

SSP/IPD 

SSP/IWD 

SSP/IDC 

SSP/GFTERM 

SSP/HYMOK 

SSP/NLVC 

SSP/PTNLVC 

SSP/MODCUP 


Mode  shape  data  print 
Particle  damping 
Wall  damping 
Distributed  combustion 
Flsndro  term 

Effect  of  high  Mach  numbers 
Nonlinear  velocity  coupling 
Nonlinear  velocity  coupling  print 
Mode  coupling 


Although  SSP/TABUL  was  deleted,  UBAR  is  printed  as  MACH  in 
INPUT  and  the  effects  of  mode  shape  may  be  examined  by  examination  of 
the  partial  integrals  DA  PL,  DAPE,  DALV  and  DAFT  added  and  printed 
in  PTL3TB.  Particle  damping  and  dietributed  combustion  were  deleted 
on  the  arguments  that  most  current  tactical  rocket  motors  employ 
reduced  or  minimum  smoke  propellants  with  small  concentrations  of 
particulates.  Wall  damping  was  also  deleted  for  small  effect:  most  tactical 
motors  have  minimal  non- burning  surfaces.  Mode  coupling,  high- Mach 
effects  and  the  Flandro  correction  were  dropped  because  of  their  contro¬ 
versial  status.  Nonlinear  velocity  coupling  was  deleted  on  the  basis 
that  the  nonlinear  effects  would  be  negligible  for  motors  with  optimised 
(high)  stability.  However,  in  connection  with  this  last  decision,  the 
linear  coupling  (now  calculated  in  STBINT)  was  changed  to  calculate  nen- 
sero  partial  integrals  only  for  sections  with  average  velocity  greater  than 
the  erosive  burning  threshold  velocity  specified  in  the  ballistic  analysis'^. 

The  logic  flow  of  E488M2  differs  little  from  SSP/MAIN.  The 
remaining  stability  integrals  are  all  calculated  in  STBINT.  The  response 
calculations  have  been  removed  from  PTLST^  and  now  constitute  subroutine 
RSPNSE. 


The  job  stacking  capability  of  SSP  has  been  replaced  with  a  formal 
time  loop  (K  =  1,  NTIMES).  The  internal  variable  NTIMES  is  selected  by 
the  calling  program.  For  the  initial  and  final  stability  calculations  with  full 
print-out,  NTIMES  =  20.  The  times  are  internally  selected  to  furnish  the 
data  required  for  stability  analysis  (principally  via  common/TIMDAT/) 
at  five- percent  intervals  of  propellant  weight  burned,  from  0%  through  9?%. 


Erosive  burning  threshold  velocity  determined  from  critical  Mach  number 
(MCRIT),  which  is  user  input  or  internally  calculated.  MCRIT  corresponds 
to  code  internal  designation  of  KR3. 
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The  use  of  weight  burned  avoid*  the  requirement  for  pre- knowledge  of 
burn  times  end  thereby  evaluates  all  designs  on  a  more  consistent  basis. 

This  tactic  provides  thorough  mapping  of  stability  parameters  for  any 
pressure-time  history  without  requiring  user  input.  For  PATSH  usage, 

NTIMES  is  reduced  to  4,  with  the  times  selected  to  furnish  the  internal 
stability  inputs  at  weight-burned  fractions  of  5%,  25%,  55%,  and  95%.  The 
uneven  spacing  provides  even  coverage  for  either  level- thrust  or  boost- 
sustain  designs.  This  spacing  also  matches  the  theoretical  trends  of 
the  individual  stability  integrals  due  to  the  fractional  change  ratio  of  internal 
volume,  port  area  or  Mach  number  by  shifting  the  time  points  earlier  in  burn 
where  the  change  rates  are  greatest. 

The  usual  error-squared  type  of  penalty  function  usually  associated 
with  PATSH  has  been  replaced  with  the  one-sided  penalty  function  OBSTAB, 
shown  in  Figure  87.  The  argument  HMIN  is  calculated  in  PTLSTB  as  the 
minimum  value  of  H  for  all  modes  analysed  at  all  times  during  burn, 
where  H  a-  ALPHA/(2*PI*FREQUENCY).  In  this  context,  H  positive  is 
good  (stable)  and  H  negative  is  bad  (unstable).  Recent  work  on  flow-driven 
oscillations  uses  H  as  a  measure  of  resonant  gain  (Ref  21),  and  therefore,  of  the 
amplitudes  of  oscillations  as  well  as  of  stability.  The  penalty  function  OBSTAB 
is  divided  into  three  regions:  unstable  (HMIN  <  0.00001),  extremely  stable 
(HMIN  >0.2)  and  the  region  in  between  where  optimisation  is  useful.  The 
unstable  region  is  clearly  unacceptable  and  large  values  of  OBSTAB  will  drive 
the  design  rapidly  toward  stability. 

In  the  extremely  stable  region,  improved  stability  would  only  require 
compromise  of  the  other  parameters,  so  OBSTAB  remains  constant  at  sero. 

A  hyperbolic  function  in  between  rapidly  drives  the  design  to  more  than 
neutral  stability,  but  with  rapidly  decreasing  stress  on  other  parameters 
as  stability  approaches  "rock  stable". 

IN  PUT  (K) 

The  input  subroutine  performs  five  functions 

o  Read/write  user  inputs  via  namelist  STABIN 
o  Default  and/or  diagnose  inconsistent  STABIN  inputs 
o  Transcribe  the  internally  transfer red/TIMDAT/data  to 
SSP  data  arrays  without  the  time  subscript 
o  Perturb  local  pressure  and  Mach  number  to  calculate  pressure 
and  velocity  exponents 

o  Write  the  input  data  set  used  for  each  time  point 

The  read  and  write  functions  are  controlled  by  IREAD  and  IPRINT 
transferred  from  the  calling  program.  Namelist  STABIN  permits  the  user 
to  specify  the  number  of  modes  to  be  analysed  (NMODE)  and  to  specify 
response  functions.  The  default  value  of  NMODE  is  four  .  This  value  is 


213 


100,000  r; 


10000 


-0.2  -0.1 


0  0.1  0.2 

h  :n 


HMIN  *  0.  2 
l  O' 6  <HMIN<ft  2 


HMIN<10 


OBSTAB  =  0.  0 

OBSTAB  =  ~-2-r  -1.0 
HMIN 

OBSTAB  =  19999(HMIN  +  1.  00001) 


HMIN  =  Minimum  H  at  any  time  for  any  mode  during  simulation 
-ALPHA 

H  =  (2)(  n  )(FREQ) 

ALPHA  =  Sum  of  gains  and  losses 
FREQ  =  Mode  frequency  at  given  time 


Figure  87.  Stability  Penalty  Function  (OBSTAB) 
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consistent  with  the  general  observation  that  longitudinal  stability  problems, 
if  any,  will  usually  arise  only  in  the  lower  modes.  The  value  also  exceeds 
the  longitudinal  stability  analysis  requirements  usually  found  in  AFRPL 
procurements.  Although  occasional  higher  mode  oscillation  problems  have 
been  encountered,  it  is  expected  that  NMODE  =  4  will  be  adequate  for  all  but 
extreme  cases. 

The  response  function  inputs  and  defaults  are  discussed  in  sub¬ 
routine  RSPNSE.  The  pressure  and  velocity  exponents  are  discussed  in 
the  section  on  subroutine  RSPNSE. 

Although  the  formats  have  been  changed  and  additional  data  printed, 
the  printout  from  INPUT  closely  follows  SSP. 

UB  CAL  C 


This  subroutine  calculates  coefficients  of  burning  perimeter  and 
port  area  variation  with  length  and  calculates  the  volume  current  (CAPQ) 
at  the  upstream  end  of  each  element.  Except  for  deleting  the  wall  and  particle 
damping  sections  and  streamlining,  the  subroutine  is  virtually  identical  to 
SSP/ UB  CALC. 

Usage  of  SSP  for  submerged  nozzles,  with  NSEC  =  20  (and  NNOZ 
<20)  is  believed  to  have  resulted  in  numerical  error  due  to  the  (J  +  l)  subscripts 
in  the  'DO  80--  'loop  in  UBCALC.  That  problem  has  hot  been  addressed 
in  the  present  program,  in  that  submerged  nozzles  are  not  available  in  SPOC. 

MODCLC,  CROSSl ,  CRQSS2,  ACOUST 

These  subroutines  solve  the  acoustic  flow  equations  for  the  standing- 
wave  mode  frequencies  using  an  iterative  search  method.  They  are  derived 
from  SSP/ MODCLC,  SSP/CRSVDX,  SSP/CRSVD2  and  SSP/ACOUST, 
respectively.  Particle  damping  and  the  SSP  print  option  NDUMP1  were 
deleted.  The  call  arguments  were  changed  by  adding  commons.  This  was 
the  only  section  of  the  stability  analysis  in  which  appreciable  difficulty 
arose  when  single  precision  arithmetic  (2 4- bit  floating-point  mantissa)  was 
used.  Occasionally,  the  search  would  locate  a  mode  frequency,  but  be 
unable  to  satisfy  the  nozzle  boundary  condition  of  U=0  and  finally  skip 
on  to  the  next  inode.  The  actual  test  argument  (UAZERO)  was  1  x  10"^. 

When  the  te  t  was  relaxed  to  1  x  10"^,  no  more  difficulties  were  encountered. 
However,  errors  in  the  fifth  significant  figure  of  stability  integrals  were 
noted  in  changing  from  56-bit  to  24- bit  mantissas.  Consequently,  the  user 
of  computers  with  more  than  24- bit  mantissas  should  examine  the  effect  of 
restoring  UAZERO  to  its  SSP  value.  UAZERO  is  specified  in  a  data  statement 
at  the  beginning  of  CROSS2  and  is  accompanied  by  appropriate  comments. 


INTGL.  G2F2FG 


These  subroutines  calculate  the  mode  amplitude  integrals  required 
in  STBINT  to  calculate  the  stability  integrals.  INTGL  was  derived  from 
SSP/INTGRL  by  deletion  of  wall  and  particle  damping  integrals  and  calculating 
several  SSP  function  subprograms,  in  the  code  stream.  The  remaining 
SSP  function  subprograms  were  combined  into  a  single  subr  outine  G2F2FG. 

STBINT,  SB  PH  AT,  SBQ  BAR 

Subroutine  STBINT  combines  the  SSP  stability  integral  subroutines 
SSP/TPCLAT,  SSP/IPCEND,  SSP/ILNVC,  SSP/IFT  and  SSP/IND.  Sub¬ 
routines  SBPHAT  and  SBQBAR  are  very  mild  revisions  of  SSP/SBPKAT 
and  SSP/SBUBARr  r ebpectively.  As  a  group,  these  subroutines  calculate 
and  store  each  required  stability  integral  over  the  length  of  each  geometry 
section.  This  tactic  permits  variation  of  the  combustion  response  functions 
along  the  length  of  the  motor  (in  subroutine  RSPNSE)  and  also  permits 
examination  of  the  contribution  of  each  section  to  the  overall  stability  margin. 
The  partial  integrals  are  printed  in  subroutine  PTJLSTB. 

RSPNSE.  ABMODL,  INT4D 

Subroutine  RSPNSE  was  created  to  calculate  combustion  response 
in  each  section  of  the  motor.  It  replaces  the  tabular  input/ interpolation 
method  included  in  SSP/PTLSTB.  Four  combustion  response  models  may 
be  specified  by  the  response  option  IRSPNS  in  namelist  STABIN: 

IRSPNb  =  1  Tabular  input  as  in  SSP 

2  Analytical  model  due  to  Culick  (Ref  22) 

3  Empirical  model  due  to  Cohen  (Ref  23) 

4  Empirical  model  due  to  Hessler  (Ref  24) 

Provisions  have  been  made  for  insertion  of  a  user-defined  response  model; 
instructions  are  included  in  comments  at  the  end  of  RSPNSE.  The  user’s 
model  \TOuld  be  inserted  into  RSPNSE  and  the  permissible  IRSPNS  t^st 
changed  in  INPUT. 

The  tabular  input  option  assumes  that  response  is  constant  in  all 
sections  of  the  motor  with  non-zero  RBAR.  The  tabular  input  in  namelist 
STABIN  is  used  to  interpolate  values  of  RPLAT,  B  PEND  and  RV  at  the 
ir  -J,5  frequency  with  FRES  as  abscissa.  This  was  the  only  option  in  SSP. 

1  interpolation  subroutine  SSP/INT4D  was  replaced  with  ITERP1  to 
provide  extrapolation.  This  was  necessary  to  avoid  error  when  the  user- 
supplied  FRES  did  not  span  all  calculated  mode  frequencies. 

Response  options  2,  3,  and  4  assume  that  the  normalized  combustion 
response  (CRPOVN)  varies  as  a  function  of  burning  rate  (RBAR)  along  the 
length  of  the  grain.  Although  only  pressure- coupled  response  was  modeled 
in  References  22  and  23,  the  models  were  extended  to  velocity- coupled 
response  using  the  methods  of  Reference  24. 


The  essential  assumption  for  the  extension  is  that  the  combustion 
response  to  heat- flux  variations  is  largely  controlled  by  the  solid  propellant. 
Representing  the  combustion  response  as  a  complex  variable,  it  follows  that 
the  response  normalised  by  the  appropriate  exponent  is  constant,  regardless 
of  the  source  of  the  heat- flux  variations.  Reference  24  proposed  that  the 
appropriate  exponents  satisfy  the  equation. 

n  n 

r  =  a(P)  P  (U)  U 


and  could  consequently  be  evaluated  from  erosive  burning  expressions. 

This  notion  was  implemented  in  INPUT  by  perturbation  of  the  burning 
rate  subroutine  (RATESB  or  USERRB)  used  in  the  ballistic  analysis 
for  the  pres  sure /velocity  field  in  each  section  of  the  motor.  The  resulting 
exponents  were  added  to  the  printout  in  PTLSTB  as  ENP,  ENU  and  (for 
end-burning  surfaces)  ENO. 

It  should  be  not*  a  that  the  perturbation  necessarily  defines  the 
existence  and  magnitude  oi  velocity  coupling  threshold  (VT)  if  the  burning 
rate  subroutine  include*?  ■»  thteshold  (either  explicit  or  implicit)  for  erosive 
burning.  Subroutine  jRA  IE..W,  #or  example,  contains  two  thresholds.  One 
threshold  is  explicit  (KT<3)  and  on  ;  implicit,  depending  on  the  relationship 
between  local  pressure  (P),  acoustic  velocity  (A),  and  the  burn  rate  inputs 
KR1,  KR2,  KR5,  and  KRbt1) 

VT  =  A*((KR1  /KR5)**(1/KR6  ))*(P**((KR2/KR6)-  1 )), 

The  lower  of  the  two  threshold  conditions  controls. 

Reference  24  points  out  that  the  velocity  response  defined  in  that  fashion 

(RU)  differs  from  the  velocity  reeponse  defined  for  the  SSP  stability  analysis 

(RV)  by  a  factor  of  the  Mach  number: 

RV  =  RU/XMACH 

Although  turbulence  effects  on  velocity  coupling  are  not  explicitly  included, 
this  method  does  include  turbulence  effects  implicit  in  the  erosive  burning 
representation  in  the  burning  rate  subroutine. 

The  final  step  in  extending  the  models  from  pressure  coupling  to 
velocity  coupling  is  use  of  the  observation  that  the  pressure- coupled  response 
function  required  by  the  SSP  analysis  mathematics  is  the  real  part  of  the 
complex  combustion  response  due  to  pressure  fluctuations.  Similarly,  the 
velocity- coupled  response  function  required  is  the  imaginary  part  of  the 


(1) 


KR1  =  A70,  AHI,  or  A.U0,  depending  on  grain  temperature;  KR2  =  XN; 
KR5  =  MPCOEF;  KR6  =  MPEXP 
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complex  combustion  response  due  to  velocity  fluctuations.  The  combined 
representation  in  RSFNSF,,  consequently  predicts  interactive  effects  between 
pressure  coupling  and  velocity  coupling  due  to  the  variation  of  burning  rate, 
pressure  and  velocity  along  the  length  of  the  motor. 

Response  option  2  uses  the  two- parameter  analytical  combustion 
response  model  form: 

Rp/n  =  AB/(  X  +  (A/X)  -  (1+A)  +  AB) 

derived  by  Culick  from  the  work  of  several  modelers  (Ref  22).  Lambda  (X) 
is  a  complex  variable  function  of  the  dimensionless  frequency  (OMEGA,  O) 
determined  by  frequency  (F),  thermal  diffusivity  (DFUSVT)  and  burning 
rate  (RBAR): 

OMEGA  =  2.  *  PI  *  F  *  DFUSVT /(RBAR  *  REAR). 

The  model  is  implemented  in  subroutine  ABMODL  using  complex  variables. 
For  user  computers  with  Fortran  compilers  that  may  not  support  complex 
arithmetic,  the  needed  real- variable  coding  is  also  supplied  in  comment 
cards  in  ABMODL  with  comment  notes  for  implementation.  Response  option  2 
requires  specification  of  the  parameters  (APARAM  and  BPARAM)  and 
DFUSVT  in  namelist  ST  A  BIN. 

Response  options  3  and  4  each  assume  that  the  magnitude  and  frequency 
of  peak  response  are  related  to  ammonium  perchlorate  (AP)  oxidizer  sizes 
(DIAAP),  tal  AP  concentration  (CONCAP)  and  concentration  of  each  size 
(CONCD): 

MPEAK  =  FACPK  *  CONCAP  *  DIAAP  **  PWR  OF  D 
FPEAK  =  F CONST  *  RBAR /DIAAP. 

It  is  also  assumed  that  the  relation  between  F  and  OMEGA  is  constant,  so 
OMEGA  =  OMEGPK  *  F/FPEAK. 

The  summation  method  used  in  Reference  24  for  multiple  oxidizer  modes 
has  been  used  for  both  models  to  assure  proper  response  at  both  zero  and 
infinite  frequency.  As  both  models  were  based  on  prtssure-coupled  T-burner 
data,  the  real  part  of  the  response  (RPOVNR)  is  used  for  scaling,  for  example, 

ARPLI  =  ARPLI  +  TEMPR  *  ENP,  where 

TEMPR  =  CONCD  *  MPEAK  **  XY,  and  the  exponent 

XY  =  FACXY  *  In  (RPOVNR) 

Reference  24  arbitrarily  set  APARAM  =  14,  and  selected  BPARAM  to  fit  the 
required  curve  shape.  Using  this  method,  two  values  of  BPARAM  were 
required  to  fit  the  curve  shape  of  Reference  23  Th«  resulting  empirical 
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)  constants  for  the  two  me  dels  are  listed  in  Table  8, 

PTLSTB  (J,  K) 

This  subroutine  performs  the  same  functions  as  SSP/ PTLSTB.  The 
code  was  rewritten  deleting  particle  damping,  wall  damping,  distributed 
combustion  and  the  Flandro  term.  In  its  present  form,  PTLSTB  sums  '.>e 
products  of  combustion  response  and  the  appropriate  oartial  integral  along 
the  length  of  the  motor,  sums  all  calculated  gains  and  losses  to  determine 
the  total  alpha  (AL),  calculates  the  fraction  of  critical  damping  (H)  and  HMIN 
and  writes  the  output  for  that  mode. 


TABLE  8 


EMPIRICAL  CONSTANTS  IN  COMBUSTION  RESPONSE  MODELS 


BASED  ON  AP  CONTENT 


J* 


Cohen  Model  Hessler  Model 


IRSPNS  = 

3 

4 

APARAM 

14. 

14. 

FCONST 

6.  0 

2.25 

PWROFD 

-1.0 

-0.  1 

FACPK 

(F  <  FPEAK) 

1056.  0 

(F  >  FPEAK) 

9.5 

BPARAM 

0.9 

1.  3 

1.4 

OMEGPK 

12.  069 

12.  621 

12.  706 

FACXY 

0.  764676 

1. 42049 

1.  57824 
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WEIGHTS 


All  weights  are  calculated  directly  from  dimensions  that  are  user 
input  or  are  internally  generated,  except  for  the  following  (which  do  not  vary 
during  a  given  problem); 

WOTHER  (components  credited  to  the  propulsion  unit) 


o 

Environmental  closure 

o  Wing  lugs 

o 

Safe-and-arm  device 

o  Igniter 

o 

Launch  lugs 

o  Miscellaneous 

WNP  (missile  components  not  credited  to  propulsion  unit) 

o  Warhead  o  Wings 

o  Guidance  and  Control 

These  two  weights  must  be  determined  by  the  user  and  supplied  as  fixed 
inputs. 


Thrust  skirt  weights  also  do  not  change  during  a  given  problem,  but 
the  forward  skirt  length  is  needed  to  calculate  overall  motor  length,  and  so 
the  code  will  calculate  forward  and  aft  skirt  weights  from  input  lengths  and 
thicknesses.  An  alternate  approach  is  to  include  Bkirt  weights  as  part  of 
WOTHER,  input  skirt  thicknesses  as  zero  and  input  an  appropriate  value 
for  forward  skirt  length;  thus  the  calculated  skirt  weight  will  be  zero  but  the 
length  is  available  for  motor  length  calculations. 

Provisions  are  made  to  include  weight  allowances  for  joint  flanges, 
even  though  there  is  no  direct  joint  design  and  weight  estimate  made  in  the 
code.  For  the  pressure  vessel,  this  allowance  is  achieved  by  adhering  to 
the  designer's  rule-of-thumb  that  the  volume  of  metal  removed  to  provide  an 
opening  in  the  closure  should  be  replaced  in  the  flange  surrounding  that  open¬ 
ing  for  attachment  of  the  mating  component.  Thus  for  the  aft  ellipse  idal 
closure  (Aft  Closure  Type  1),  the  weight  of  the  center  segment  shown  in 
Figure  88  is  allocated  to  the  aft  flange  on  the  case  and  is  separately  identified 
as  WFLGA.  An  allowance  for  the  weight  of  the  flange  on  the  nozzle  is  pro¬ 
vided  by  virtue  of  the  nozzle  computation  scheme  that  measures  all  nozzle 
thicknesses  normal  to  the  interior  surface.  Therefore,  there  is  a  portion  of 
the  nozzle  that  can  rightly  be  designated  as  "flange"  (Figure  88  )#  although 
it  is  not  identified  as  such  in  the  code  output.  At  the  forward  end,  the  closure 
weights  are  calculated  as  if  there  were  no  opening,  thus  preserving  the  above 
quoted  rule-of-thumb,  bv*$  as  with  the  nozzle,  there  is  no  special  identification 
of  the  forward  flange  weights. 

There  is  a  small  volume  of  material  accounted  for  twice,  once  in  the 
skirt  volume  and  once  in  the  ellipsoidal  closure  volume  (Figure  89  ).  There 
is  always  a  fillet  between  the  skirt  and  closure  and  the  duplicated  volume  makes 
at  least  some  allowance  for  this  extra  material. 
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Allowance  For  Case 
Aft  Flange  Weight 


Figure  88.  Weight  Considerations  at  Ait  Case  Opening 
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MISCELLANEOUS 
Liner  and  Insulation 

Liner  ia  described  to  SPOC  aa  a  single  material  having  a  constant 
thickness  everywhere  in  the  motor. 

Two  classes  of  insulation  may  be  defined  in  the  code;  one  ia  allocated 
specifically  to  protect  the  grain  against  the  effects  of  external  aerodynamic 
heating  and  the  other  is  to  protest  the  case  against  internal  heating  by  the 
combustion  products. 

Aerodynamic  heating  insulation  is  specified  through  input  of  a 
thickness  (TAEROJ.)  and  density  (DELAI).  This  insulation  has  a  constant 
thickness  over  the  entire  cylindrical  portion  of  the  case;  there  are  no  pro¬ 
visions  for  aerodynamic  heating  insulation  in  the  closure?.  The  thickness 
cannot  be  adjusted  during  a  given  machine  submission. 

Internal  heating  insulation  is  specified  through  input  of  a  number  of 
thicknesses  and  a  single  value  of  density  (DELINS).  This  insulation  haB  a 
constant  thickness  around  the  periphery  at  any  given  longitudinal  station, 
but  the  thickness  can  vary  along  the  length  of  the  motor  in  accordance  with 
the  input  feature?  provided  for  the  different  grain  types.  In  particular, 
this  means  that  he  internal  insulation  cannot  be  contoured  to  match  a  star 
or  wagon  wheel  grain;  instead,  the  user  must  select  an  average  thickness 
with  which  to  account  for  the  volume  displaced  by  the  insulation.  Internal 
insulation  thickness  is  not  adjusted  during  a  given  machine  submission;  the 
user  must  make  a  selection  for  the  initial  inputs  at  the  different  locations 
in  the  motor  based  on  his  knowledge  of  what  the  heating  conditions  are 
expected  to  be.  After  reaching  a  full  or  partial  solution  to  the  problem,  these 
thickness  choices  are  reviewed,  based  on  the  latest  definition  of  heating 
conditions,  and  adjustments  to  the  thicknesses  are  made  if  deemed  appro- 
pi’ ate  (or  if  dictated  by  thermal  analyses  made  external  to  the  optimization 
problem);  this  is  the  basic  preliminary  design  approach. 

Insulation  in  the  Type  1  (ellipsoidal)  forward  and  aft  cloaures  is 
described  with  a  minimum  thickness  at  the  case-closure  tangent  points  and 
a  maximum  thickness  at  the  case  opennings.  Insulation  on  the  flat  plate 
forward  closure  Type  2  os  Type  3  has  a  separate  thickness  input.  All 
internal  heating  insulation  has  the  sa  ne  density;  a  different  density  may 
be  input  for  the  aerodynamic  heating  insulation. 

Stress  - R elief  Boots 


Stress  relief  boots  may  be  specified  through  input  of  a  thickness 
at  the  case  openings  and  a  separate  density  value.  The  boot  tapers  to  zero 
thickness  at  the  case-closure  tangent  points.  Stress  relief  boots  (or  flaps) 
cannot  be  defined  for  closure  Types  2  or  3,  but  allowances  can  be  made  for 
them  by  adjusting  the  insulation  thicknesses. 
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USER  MODELS 


Provision*  have  been  made  for  the  uaer  to  supply  special  mathe¬ 
matical  models  to  describe  certain  characteristics  of  the  motor.  They  are: 

(1)  Propellant  burn  rate 

(2)  Propellant  nominal  strain  endurance 

(3)  Propellant  rheological  property 

(4)  Motor  costs 

(5)  Impulse  efficiency 

(6)  Combustion  response 

The  procedure  for  employing  a  us*r- supplied  model  is  as  follows: 

(1)  In  the  subroutines  that  are  furnished  (Table  9.  ) 

(a)  Code  the  mathematical  model,  making  sure  the  dependent 
variable  has  the  nomenclature  given  in  Table  9. 

(b)  Add  any  common  statements  required  to  furnish  the  inde¬ 
pendent  variables  required  by  the  model. 

(c)  Furnish  any  WRITE  commands  that  are  needed  to  print 
information  from  the  model.  See  further  discussion  below. 

(2)  Compile  the  subroutine  and  link  with  the  remainder  of  the  code. 

(3)  Set  the  flag  in  the  namelists  shown  in  Table  9  to  show  that  a 
given  user  model  is  being  furnished. 

The  internal  flag  of  IPRINT  is  used  to  control  when  computation 
results  are  printed. 

IPRINT  ~  1  All  write  statements  are  executed 

IPRINT  =  0  Only  PATSH  output  is  printed 

Therefore,  WRITE  commands  in  the  user-model  subroutines  should  be 
structured  such  that  they  are  executed  when  IPRINT  =  1  and  are  branched 
around  when  IPRINT  =  0. 

If  data  must  be  supplied  to  the  model  by  the  user,  a  READ  command 
must  be  included  for  a  namelist  defined  by  the  user.  Common  TRIGR 
contains  the  flag  IPRINT  (discussed  above)  and  a  similar  flag  to  read  data 
(IRE  AD). 

IREAD  =  1  To  read  data 

IREAD  =  0  To  bypass  read  command 

IREAD  =  1  only  for  first  pass  through  COMP. 

IPRINT  =  1  for  both  the  first  and  last  passes  through  COMP. 
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TABLE  9 


PROPELLANT  BURNING  RATE 


The  user- supplied  propellant  burn  rate  model  muet  be  installed  in 
subroutine  USERRB  The  flag  RBMDL  *  T  must  be  input  in  nameliet  CCMTRL 
so  that  USERRB  will  be  called  whenever  burn  rate  muet  be  calculated.  The 
common  statement  /IB/  must  be  employed  In  subroutine  USERRB.  to  furnish 
values  for  pressure  (P),  Mach  number  (MACH),  and  rate  scale  factor  (SF) 
and  to  pass  back  the  calculated  burn  rate  (RATE).  If  the  user  wishes  to 
identify  which  burn  rate  equation  out  of  several  possibilities  is  being  employed 
in  subroutine  USERRB,  the  common  statement  /SEC3/  must  be  included  to 
furnish  the  variable  RATEQU;  values  for  RATEQU  will  be  printed  out  at 
time  zero  in  the  ballistic  simulation  for  each  plane  along  the  grain. 

The  burn  rate  subroutine  is  called  from  subroutines  HITEMP, 

LOTEMP  and  ONETMP  with  a  pressure  averaged  over  the  ballistic  simulations 
at  high  temperature,  low  temperature  and  single  temperature,  respectively. 
The  purpose  of  the  calculations  at  this  point  in  the  code  is  to  obtain  an 
average  burn  rate  without  cross  Row  effects  to  send  to  the  SPP  impulse 
efficiency  model.  Therefore,  Mach  number  is  not  supplied  from  these 
subroutines,  and  a  flag  is  set  (IFLAG  =  1)  to  mark  where  the  call  to  the 
rate  subroutine  originated.  Then,  in  the  rate  subroutine,  IFLAG  =  1 
triggers  a  return  to  the  calling  routine  after  a  non-erosive  rate  has  been 
calculated.  The  user-supplied  subroutine  USERRB  must  contain  this  same 
response. 

Other  commons  that  already  exist  may  be  included  in  USERRB  in 
order  to  furnish  the  data  needed  for  the  user  model. 

Temperature  effects  on  burn  rate  can  be  included  through  the  use 
of  TMPUR  as  temperature  in  the  neW  model.  The  variable  TMPUR,  which 
is  contained  in  the  common  statement  MODELS,  is  equated  to  THI  just  prior 
to  the  high  temperature  ballistic  simulation;  it  is  then  equated  to  TLO  just 
before  the  low  temperature  simulation.  Thus,  TMPUR  has  the  appropriate 
value  at  the  time  USERRB  is  called  in  the  ballistic  simulation. 

A  typical  USERRB  is  given  in  Table  10. 

PROPELLANT  STRAIN  ENDURANCE 

The  user-supplied  propellant  strain  endurance  model  must  be  in¬ 
stalled  in  subroutine  USERSE.  The  flag  SEMDL  =  T  must  be  input  in  namelist 
STPROP  so  that  USERSE  will  be  called  out  of  subroutine  PROPST.  The 
model  must  furnish  the  nominal  strain  endurance  (SENOM)  for  comparison 
with  thermally  induced  strains.  The  nominal  value  is  devalued  for  mix-to- 
mix  variations  and  aging  degradation  in  PROPST.  SENOM  is  returned  from 
USERSE  through  the  catling  argument. 

If  the  model  is  not  furnished,  SEMON  is  a  constant  user-supplied 

input. 
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TABLE  10 

TYPICAL  USER-SUPPLIED 
BURN  RATE  SUBROUTINE  (USERRB) 


ISN  000? 
ISN  0004 

ISN  000*. 


SUBROUTINE  ustnco 

COMMON/TS  IC.R/IRE  AO,  IPF  |M  ,10111  I ,  I  QUIT  I ,  ir  LA.', 
COMMON  /r»t  I  ST  /  C(»S|,ASMTCHIT';i,nAJ,MlNtAT,I.TOTAl 


COMOOI?.! 

rnMnni*.n 

crMnniin 


COMMON  /tB/APHIB.  .SPNIO,  SFHIB,  SFNI6*  COAHIB,  CO*  N  IB,  PI  ASIA,  fOMOn’MI 
.  ®HIA,  TH1A,  OITMIB,  (MIS,  WPTM1A,  RTEM1H,  P,  T,  DTI  T  A,  CONOOINO 

.  Ui  wnnr,  mach,  sf,  coaatf,  bate,  scNcnr  coaoo>f.o 

comootto 

COMMON/ NODEIS/TM®!)* 

COMMnN/BFQlR/rpT.PNNAX.MACMAX,  AE,PATNI  ,CVT  ,r.VTH,fVl  ,r  VP  I,  r  MIN,  COMOOI  10 
.  I TACQ.TANXLM,  TANNIN,  MAfUM, 4  A,  P IGN  ,CV»,FSUU  ,  T  Stl.P ,  FM.AX  .NN**,  ACl  IM.CONOOIAn 
.  PMXl  IN,  put  t  ,PY  I  El  0  ,T  1C ,  MCOCAS*  PAC  AS,OEl  VAQ,11TAQ»VTRQ»V80,ISPT0,  CPH-IO  ISO 
.ISPMI .ISPIO.NTFAPS.ACI  PAX,NrTNEA,TH| ,1TMIN  rnMOO'M) 

roMOoifn 


isn  on in 


ISN  0011 
ISN  001? 


ISN  001 » 
ISN  001  A 
ISN  OOIS 
ISN  0016 


r  CNNilN/APO  AT  A/CGNf.AP  ,  l  •'TOE  ,P  1  A  API  1 1  ,  CllNCH  I  SI 

COMMON  /CMEMIN/OINO.FUn  ,OXA  II  l,OXA  IS  I.WCAU,  ACATS  .  SI  AB,  I  01  Nn, 
.IFUtl.lOXAt II,  ICXBOI.I  ’’CL,  tACS.ISCl,  I STAB.FUEl AX , PCI  MX . AC SN* , 
•  TSMAX.UA.LIL 

COMMON  /SEC  4/orif- ,  GAMMA,  BOAS,  TS,  RA1,  AH?,  API,  AAA,  APS, 
.AAA.USONIf ,CIA|,CIB?,C I  SI, AA IF  00 

®F  At  MACH, API  AP?,ARS,1PA,APS,AR6 
AFAl  NC,inflAAl,lOnAAJ,lPr,ATO,LPBTO 

USF"  MUST  INPUT  VAIUFS  FCR  MA*  ANO  MIN  BM»N  » ATF  CONSTRAINTS, 
TEMPEP.ATUAE  COEFFICIENT  OF  CSTAS  (MCI, 

TFNPFRATUAE  COEFFICIENT  OF  PRESSURE  IPIKI, 

BEFORE  CONAIUNG  THIS  SUBROUTINE 

RBMAX  •  100. 0 
ROM  IN  •  0.0 
PIA  •  O.OOIS 

MC  •  o.ooos 


nno.vmo 
ooo no  ’on 

00000? TO 


ISM  001 7 
ISN  OOIA 
ISN  0010 
ISN  00 ’0 
ISN  00?1 
ISN  00?? 


CONVERT  INPUTS  TO  UNI  Tr  NtFOFD  FOR  »ATC  NOOrt  ICIIANCT  WT  t  CHI 
FRACTION  TO  HEIGHT  AEPrSNTI 

API  »  OXAI  II  •  100. 
ap>  •  rxA(?i  •  ino, 

API  •  OXAI II  *  100. 

TOTAP  »  ICIXAIII  ♦  OXAI  ?  1  ♦  OXAI  III  ♦  100. 

Al  •  FUEI  *  ion. 

FTC  *  PC  AT  S  •  ICO. 


CALCUt ATE  CONSTANT  TCP'S 


ISN  0021 
ISN  0026 
ISN  0025 


08 AR I  *  I AP I *01 AAP | t I  *  AP?*OI AAP I ? I  »  AP l*OIAAP|  U l/TUTAP 
LOAARI  -  AECG10(10.*0BA°|| 

LOOARI  *  I  AP1*AEOGIO|  l  U,*OI  A  API  1 1 1  ♦  AP?«  AEOCil  01 1  n  .  *01  AAPt  '  II  , 
.  API*  AinC.|OUi'.*OIAAPI  II I  l/TUTAP 

te“pi  »  o.sbt*u:aa«i 
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Table  10  (Continued) 

Typical  Uaer-Supplied  Burn  Rate 
Subroutine  (USERRB) 


ISN  00/7 
ISN  007# 
ISN  007*7 


ISN  0010 
ISN  00 M 
ISN  0)1,' 


ISN  non 

ISN  001* 
ISN  00  IS 
ISN  00)7 


ISN  00)4 


ISN  00*0 
ISN  00*/ 


ISN  00*1 


ISN  00*5 

ISN  00*7 

ISN  00*# 
ISN  00*4 
ISN  0050 
ISN  00*1 
ISN  00*2 
ISN  0051 
ISN  005S 


ISN  0057 
ISN  005# 


ISN  0054 


TfNP2*-0.00)?t*lOf  AR*lOBARi 
TERR)  »  O.OS*A»fEnUOBARl 
IE"#*  «  -0*0221 RALRFEO 

c 

C  CALCULATE  BU«N  RATC  AAAANETEBS  AT  TO  OTC.  #ATC  IN  NOOFI  IS 
C  SCALED  TO  10  TINES  TRUE  BATE 
C 

L0f.A70  ■  If  NPl  ♦  TEH#2  ♦  TtRP)  ♦  TF  #"*  -  | .  25 
ATO  »  I  10. •»  L0GA7OI/10. 

AN  •  0*01  36*701 AP  *  O.QOOl A2*f EC *TOT AR*  AL  -  0 .024R*TFI)*ff  O  - 
.  0,255*L08AR)  ♦  0.000266*AL*AL 

C 

C  CHECK  KUAN  NAIE  AT  70  DEC,  1000  *SIA  AOAINST  CONSTRAINTS 
C 

l  #070  «  l CO A  TO  *  XN*3.0 
#870  •10.**IL#r  .  |.ot 

iriKB70.GT.A8N.  <1  OARBRX  ■  I  IRB70  -  AB*AX > **2 1  * l.Ot A 
IFIRnTO.LT.RBNINI  OBABNN  •  IIABTO  -  * ONIN ) »*2 I  *  1 . OF  * 

C 

C  UPDATE  PENALTY  SUNHAT1CN 

OBJ  »  OBJ  ♦  OBHBNX  *  OBABNN 
C 

C  CALCULATE  PRESSURE  COEFFICIENT  AT  TFNPEAA TUBE  tXTRFNFS 
C 

IE(NTEN*>S.NE.2I1NPUA  «  TO. 

AiHPup<Aro«rxpiiPlx*ll.o-xNt-Hci  •  itnpua  -  to.oii 

c 

C  WAITE  BASIC  DATA  ETON  PATE  NOOEL 

C 

C  PAINT  ON  FIRST  CALL  TQ  USEAA8  DURING  FIRST  PISS  THROUGH  CONP 

IF  I  IPRINT.EQ.  1.  ANO. IREAD.EO.l. ANO. IHRST.NE.IIGQ  TO  20 
C 

C  PRINT  ON  FIRST  CALL  TO  USFAR8  DURING  LAST  TASS  THROUGH  COMP 

IFIIPAINT  .E<1. 1  .ANn.IREAO.EQ.fi.  AND.  ILAST  .NE.l  >G0  TO  20 
C 

GO  TO  30 
C 

20  WR I TE I  6  *  1000 1 

WR1TFI  A,  1001  IRnTO,  XN,  ATO.OBRBHX .OBRBNN 
WRITF <6,10021  API, APJ.AP 3, TOT AP.FEO.Al 
WR  ITE(6,lOQ3)OBARl tLORAAl, LOBAR) 

WRITE  1 6, IODAILUGA70, ATNPUA , TNPUA 

I F I  IPRINT.EQ. I,  ANT).  IRF  AO.EO.  I.  AND.  IF  IRST  .NF  .  1 1 1 F  IR  ST.  1 
in IPRINT.EQ. I, ANO. IRE  AO.EQ .0. ANO. ILAST.NE.il ILAST*l 
C 

C  CALCULATE  BURN  RATE  AT  CURRENT  CONDITIONS 
C 

30  RATE  *  SE  •  ATMPUAMP**XNI 
AATEOU  «  1.0 
C 

C  IF  1FLAC*1,PAT£  SUBROUTINE  HAS  BEEN  CALi-EO  FRIIH  SUBR0UT1NF  HIUHP, 

C  LOTENP.OR  ONETNP  TC  CALCULATE  RATE  AT  AN  AVERAGE  PRESSURE  WITHUUI 

C  EROSIVE  BURNING.  NUST  HAVE  THIS  LINE  CF  COOE  IN  ANY  RATF  SUA*0>|HNr 

C 

IFI IELAG.EG.il  RETURN 
C 

C  CDRATE-l.O  FOR  NO  FRQSIVE  BUR  INC,  -2.0  FUR  EROMVE  BURNING.  M  IS 

C  USEO  TO  OBTAIN  INTEGER  TO  TEST.  N»1  CAUSES  »r  USN  WITHOUT 
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Table  10  (Continued) 

Typical  User-Supplied  Burn  Rate 
Subroutine  (USERRB)  ~~ 


ISN  0061 
ISN  0062 


ISN  0066 
ISN  0066 


ISN  OObS 
ISN  0069 
ISN  007 1 
I  SN  0071 
ISN  0076 
ISN  0076 


ISN  0070 
ISN  0080 
ISN  0081 
ISN  0081 

ISN  0085 
ISN  0006 

ISN  0087 

ISN  0088 
ISN  0009 
ISN  0090 
(SN  0091 


C  CALCJLAt ING  EROSIVE  BURNING  CfFECTS.  MUSI  HAVE  THFSF  TWO  L INF S 
C  Of  COOE. 

c 

M.  *  CORATC  ♦  0.6 
IFIN.FO.il  UPTURN 
C 

C  KR 1*NCR  IT  ,KP4aXM*X, MACH* LOCAL  MACH  NUNRER  .  CAUSES  SADERHDLM  MOOH 

C  TO  OF  SKIRPfS  If  K»1  CESS  THAN  2FA'  OR  IHCAI.  MACH  NUMBER  LFSS  THAN 

C  NCR  IT. 

C 

IFIXR3.LE.O.O)  GO  TO  10 
If (MACH.LE.KR31  GO  TO  10 
C 

C  F  805! VE  DURING  WITH  SADERHOIM  MODEL 

C 

X  *  At OGIO. 06 768*IR**0. 74/RAT G 1**0 .49481 
IFIX.LC.O.OI  X*0.0 
IFIKR4.GT. 0.0011  X  »  KR6 
RATER  «  RATE»IMACH/><A3I**X 
IF  |PATF7.GT.RATF)RATF(}U*?.0 
I M "ATF2 .GT .R  ATF 1  FAT F*R  AT  E  2 
C 

C  FROSIVE  BURING  KITH  NACH-PRESSURE  PROOUCT  MODEL 
C 

10  IFINACM  *  P.LE.0.0 1  RETURN  -  -  -  -  - 

RATE3  *  KR5*INACH*PI**XR6 
IFIKATE1.GT .RA IE 1RATE0U*!.  0 

IT  IRATC3.GT.RATFIRATr*RATE3  -  -  -  - .  -  ■ 

C 

1000  FORMAT! IH1 ,41H  OASIC  DATA  FROM  USER  BURN  RATF  MOREL  1 

1001  FflRNATUHn,5HRB70*,El2.4,3X,lMXN*,El?.4.1XtAHA70-,F12.4, - - 

•  3X,7HOBROMX  =  ,F.  i2.4,3X,  7H08RBMN^  .El  2.4  I 

1002  FORMAT  I IH0,4HAPI-,E12.4,1X,4HAP2*,EI2.4,1X,4HAP1’,E|2.4, 
.1X,6HTOTAP=,C12.4.1X,4HFEO*,E12.4, IX, JHAL-,E12.4> 

1003  FORMAT  1 1H0, 7H0B AR l« ,C 12 .4 , 3X. THLD8 AP l« , El  2.4, IX  7HL RBAR  J* , E 1 2 . 4  I 

1004  FORMAT  1 1H0, 7HLOGA70- ,  El  2.4 , 3X,  7MATMPUR- »F  1  2.4 , 3X,  7IITMPUR*  ,Fl2.b) 
RFTURN 

END 
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PROPELLANT  RHEOLOGICAL  PROPERTY 


The  user-supplied  propellant  rheology  model  must  be  installed  in 
subroutine  USERRH.  The  flag  EOMMDL  =  T  must  be  input  in  namelist 
INGLIM  so  that  USERRH  will  be  called  out  of  subroutine  TCHEM. 

This  parameter  (EOM)  must  be  defined  by  the  user;  it  can  be  end-of- 
mix  viscosity  (hence  EOM),  shear  stress,  or  any  other  measure  of  rheo¬ 
logical  characteristics.  The  intent  is  that  this  model  be  used  when  propellant 
ingredient  concentrations  are  being  adjusted,  so  that  some  control  can  be 
exercised  on  the  propellant  processibility.  If  the  flag  has  been  set 
(EOMMDL  =  T),  the  output  of  the  model  (EOM)  is  compared  with  an  input 
maximum  limit,  and  a  penalty  (OBJEOM)  is  calculated  if  the  limit  is 
exceeded.  EOM  is  returned  from  USERRH  through  the  calling  argument. 

MOTOR  COST 


The  user-supplied  motor  cost  model  must  be  installed  in  subroutine 
USERCS.  The  flag  CSTMDL  =  T  must  be  input  in  namelist  CONTRL  no  that 
USERCS  will  be  called  out  of  subroutine  COMP.  The  model  must  furnish 
the  parameter  COST,  which  is  used  only  an  one  of  the  payoff  parameters 
The  units  can  be  either  total  project  cost  or  unit  cost.  COST  is  returned 
from  USERCS  through  the  common  statement  MISL. 

IMPULSE  EFFICIENCY 

The  user-supplied  impulse  efficiency  model  must  be  installed  in 
subroutine  USEREF.  The  flag  EFMDL  =  T  must  be  input  in  namelist 
CONTRL  so  that  USEREF  will  be  called  out  of  subroutine  COMP.  If  a  user 
efficiency  model  is  called,  then  the'SPP  model  cannot  be  called  (i.  e. , 
SPPETA  =  F  is  required).  The  model-furnished  impulse  efficiency  then  is 
used  in  the  ballistic  simulation.  Efficiency  is  returned  from  USEREF 
through  the  calling  argument. 

COMBUSTION  RESPONSE 


The  user- suppli ed  combustion  response  model  must  be  installed  in 
subroutine  RSPNSE,  which  is  the  subroutine  where  all  internal  models  are 
located.  Thus  many  sources  of  input  data  are  already  available.  Entry  to 
the  model  is  statement  number  500,  The  flag  IRSPNS  =  5  will  cause  the 
user- supplied  model  to  be  called. 

Combustion  response  is  returned  through  the  common  ALF 
already  furnished. 
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ERROR  MESSAGES 


Error  messages  are  printed  whenever  code  execution  is  terminated 
under  abnormal  conditions.  Other  messages  are  provided  the  user  to  show 
conditions  under  which  the  code  is  operating.  The  messages  that  originate 
from  subroutines  that  were  developed  for  SPOC  are  self-explanatory  and 
do  not  require  a  complete  listing  in  the  user's  manual  (except  as  discussed 
below).  Generally  these  messages  explain  what  caused  the  termination  and 
provide  data  to  show  the  abnormal  condition;  if  the  solution  is  not  self-evident, 
instructions  are  given. 

There  is  a  large  group  of  termination  messages  that  originate  from  the 
subroutines  comprising  the  ballistic  simulation  module.  None  of  these  con¬ 
ditions  should  be  encountered  because  the  purpose  of  the  SETUP  subroutine 
is  to  prevent  invalid  conditions  from  occurring  during  the  optimization 
process.  Nevertheless,  the  ballistic  module  error  messages  are  given  in 
Table  11  in  the  event  that  a  combination  of  inputs  and  optimization  searches 
produces  an  invalid  situation. 


The  messages  are  listed  alphabetically.  Those  that  begin  with  asterisks 
are  listed  alphabetically  at  the  end  of  the  table.  With  the  message  is  the 
reason  that  message  was  triggered.  The  last  column  of  the  table  is  the 
format  statement  number  and  the  subroutine  ’n  which  that  statement  is 
found.  Note  that  some  of  the  messages  speak  of  a  "Type  1"  or  "Type  2" 
grain;  these  are  designations  internal  to  the  ballistic  simulation  module 
and  have  no  direct  relation  to  the  grain  types  that  can  be  selected  by  the  user 
of  SPOC.  For  information  purposes,  the  correspondence  is 


SPOC  Grain 
1X21 _ 


Ballistic  Module 
Grain  Type _ 


1  (Star) 

2  (Wagon  Wheel) 

3  ( Finocy  1) 

4  (Conocyl) 

5  (CP) 


2 

2 

3 

1  (Fwd,  segment) 
3  (Remaining) 


Incompatible  inputs  for  the  trajectory  simulation  are  detected  and 
identified  with  the  flag  ISTGP  Definitions  of  1ST  OP  are  given  in  Table  12. 
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in  P.aliistic  Simulation  Module 


**  INVALID  ALPHA  01  *** 


Table  11  (Cont'd) 

Error  Messages  in  Ballistic  Simulation  Module 


TABLE  12 

EXITS  FROM  TRAJECTORY  SIMULATION  FOR 
INCOMPATIBLE  INPUTS  USING  ISTOP 

ISTOP _ Explanation _ 

1  Both  launch  velocity  and  launch  Mach  number  are  input 

2  ITERM  <  1,  or  ITERM  >  12 

3  F PAAL  <  -90.0,  or  FPAAL  >90.0 

4  AREF  <0.0 

5  IATMOS  <  0,  or  IATMOS  >4 

6  ALTAL  <  0.  0 

7  AEX  <0.0 

8  RNGAL  <0.0 

9  WNP  <  0.  0 

^0  Failure  to  pass  the  test  (IPRDEG  t  0  and  IPRDET  i  1) 

11  RGFPAL  <  0.  0 

12  WMI  <0.0 

13  FDELT1  <0.0 

14  FDELT2  <0.0 

16  ITERM  =  1  and  TERTIM  <0.0 

17  ITERM  =  2  and  TPHANE  <  0.  0 

18  ITERM  =  3  and  SPTERM  <  0.  0 

19  ITEP  *  '  =  4  and  RGTERM  <  0.  0 

20  ITT  '  .iv.  =  6  and  MTERM  =  0.  0 

21  ITERM  =  7  and  VLTERM  =  0.0 

22  ITERM  =  8  and  FPAAL  <-90.0 

23  ITERN  :  and  FPAAL  <  FPATRM 

24  ITERM  =  8  and  FPAAL  >90.  0 

25  ITERM  =  9  and  ACTERM  >0.  0 

26  ITERM  =  11  and  RGFPTM  <0.0 

27  ALTERM  <0.0,  and  ALTAL  <|aLTERm|,  and  FPAAL  >0.0 

28  ALTERM  <  0.  0,  and  ALTAL  <  ALTERM,  and  FPAAL  £  0.  0 
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